Fennec.h
Go to the documentation of this file.
1 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////
2 // Make CSV
3 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////
4 void MakeCSV(TFile* inFileData, TFile* inFileMC, std::string detector, std::string sample, std::string period, std::string outputName,
5  bool byDB, int runsPerDiv, std::string calibMethod, bool mcRun, std::string user_dir, bool verbosity){
6  std::cout << "Determining calibration...";
7 
8  // Determine number of run divisions. 0 means calibrate monolithically.
9  long timeDiv=0;
10  long loRun(0),hiRun(0);
11  int fileRunsPerDiv(0);
12  TH1D *metadataPlot;
13  inFileData->GetObject("metadataPlot",metadataPlot);
14  if( metadataPlot ){
15  fileRunsPerDiv = metadataPlot->GetBinContent(2);
16  timeDiv = metadataPlot->GetBinContent(5);
17  loRun = metadataPlot->GetBinContent(3);
18  hiRun = metadataPlot->GetBinContent(4);
19  }
20 
21  if( fileRunsPerDiv!=runsPerDiv ){
22  std::cout << "WARNING!!!!!"
23  << "\n\t runsPerDiv requested (" << runsPerDiv << ") doesn't match what file produced with (" << fileRunsPerDiv << ") ..."
24  << "\n\t Defaulting to calibrate full file at once (runsPerDiv=0)..." << std::endl;
25  runsPerDiv=0;
26  timeDiv=0;
27  }
28 
29  int tv = 1;
30  unsigned int nChan = 14;
31  if( detector.find("fd")!=std::string::npos ) nChan = 14;
32  else if( detector.find("nd")!=std::string::npos ) nChan = 4;
33 
34  std::vector<TH1D*> hmeuReco_x;
35  hmeuReco_x.resize((nChan+1)*(timeDiv+1));
36  std::vector<TH1D*> hmeuReco_y;
37  hmeuReco_y.resize((nChan+1)*(timeDiv+1));
38  std::vector<TH1D*> hmeuTrue_x;
39  hmeuTrue_x.resize((nChan+1)*(timeDiv+1));
40  std::vector<TH1D*> hmeuTrue_y;
41  hmeuTrue_y.resize((nChan+1)*(timeDiv+1));
42 
43  std::vector<double> meu_x;
44  meu_x.resize((nChan+1)*(timeDiv+1));
45  std::vector<double> meu_y;
46  meu_y.resize((nChan+1)*(timeDiv+1));
47  std::vector<double> meu;
48  meu.resize((nChan+1)*(timeDiv+1));
49 
50  std::vector<double> meuRMS_x;
51  meuRMS_x.resize((nChan+1)*(timeDiv+1));
52  std::vector<double> meuRMS_y;
53  meuRMS_y.resize((nChan+1)*(timeDiv+1));
54  std::vector<double> meuRMS;
55  meuRMS.resize((nChan+1)*(timeDiv+1));
56 
57  std::vector<double> meuErr_x;
58  meuErr_x.resize((nChan+1)*(timeDiv+1));
59  std::vector<double> meuErr_y;
60  meuErr_y.resize((nChan+1)*(timeDiv+1));
61  std::vector<double> meuErr;
62  meuErr.resize((nChan+1)*(timeDiv+1));
63 
64  std::vector<double> mev_x;
65  mev_x.resize((nChan+1)*(timeDiv+1));
66  std::vector<double> mev_y;
67  mev_y.resize((nChan+1)*(timeDiv+1));
68  std::vector<double> mev;
69  mev.resize((nChan+1)*(timeDiv+1));
70 
71  std::vector<double> mevRMS_x;
72  mevRMS_x.resize((nChan+1)*(timeDiv+1));
73  std::vector<double> mevRMS_y;
74  mevRMS_y.resize((nChan+1)*(timeDiv+1));
75  std::vector<double> mevRMS;
76  mevRMS.resize((nChan+1)*(timeDiv+1));
77 
78  std::vector<double> mevErr_x;
79  mevErr_x.resize((nChan+1)*(timeDiv+1));
80  std::vector<double> mevErr_y;
81  mevErr_y.resize((nChan+1)*(timeDiv+1));
82  std::vector<double> mevErr;
83  mevErr.resize((nChan+1)*(timeDiv+1));
84 
85  for( int iTime=0; iTime<=timeDiv; ++iTime ){
86  if( runsPerDiv>0 && iTime==timeDiv ) break; // include iTime==0 when that's the only one, otherwise break when iTime==timeDiv
87  std::cout << "Loading plots for timeDiv " << iTime << "...";
88  if( timeDiv==0 ){
89  inFileData->GetObject("th1d_nhits_resp_flatW_xview",hmeuReco_x[(nChan+1)*iTime]);
90  inFileMC->GetObject("th1d_nhits_truededx_flatW_xview",hmeuTrue_x[(nChan+1)*iTime]);
91  inFileData->GetObject("th1d_nhits_resp_flatW_yview",hmeuReco_y[(nChan+1)*iTime]);
92  inFileMC->GetObject("th1d_nhits_truededx_flatW_yview",hmeuTrue_y[(nChan+1)*iTime]);
93  }
94  else{
95  if( mcRun ) inFileData->GetObject("th1d_nhits_resp_flatW_xview",hmeuReco_x[(nChan+1)*iTime]);
96  else inFileData->GetObject(TString::Format("th1d_range%u_nhits_resp_flatW_xview",iTime),hmeuReco_x[(nChan+1)*iTime]);
97  inFileMC->GetObject("th1d_nhits_truededx_flatW_xview",hmeuTrue_x[(nChan+1)*iTime]);
98  if( mcRun ) inFileData->GetObject("th1d_nhits_resp_flatW_yview",hmeuReco_y[(nChan+1)*iTime]);
99  else inFileData->GetObject(TString::Format("th1d_range%u_nhits_resp_flatW_yview",iTime),hmeuReco_y[(nChan+1)*iTime]);
100  inFileMC->GetObject("th1d_nhits_truededx_flatW_yview",hmeuTrue_y[(nChan+1)*iTime]);
101  }
102  if( byDB ){
103  for( unsigned int iChan=1; iChan<=nChan; ++iChan ){
104  inFileMC->GetObject(TString::Format("th1d_db%u_nhits_truededx_flatW_xview",iChan),hmeuTrue_x[(nChan+1)*iTime + iChan]);
105  inFileMC->GetObject(TString::Format("th1d_db%u_nhits_truededx_flatW_yview",iChan),hmeuTrue_y[(nChan+1)*iTime + iChan]);
106  if( timeDiv==0 ){
107  inFileData->GetObject(TString::Format("th1d_db%u_nhits_resp_flatW_xview",iChan),hmeuReco_x[(nChan+1)*iTime + iChan]);
108  inFileData->GetObject(TString::Format("th1d_db%u_nhits_resp_flatW_yview",iChan),hmeuReco_y[(nChan+1)*iTime + iChan]);
109  }
110  else if( mcRun ){
111  inFileData->GetObject(TString::Format("th1d_db%u_nhits_resp_flatW_xview",iChan),hmeuReco_x[(nChan+1)*iTime + iChan]);
112  inFileData->GetObject(TString::Format("th1d_db%u_nhits_resp_flatW_yview",iChan),hmeuReco_y[(nChan+1)*iTime + iChan]);
113  }
114  else{
115  inFileData->GetObject(TString::Format("th1d_range%u_db%u_nhits_resp_flatW_xview",iTime,iChan),hmeuReco_x[(nChan+1)*iTime + iChan]);
116  inFileData->GetObject(TString::Format("th1d_range%u_db%u_nhits_resp_flatW_yview",iTime,iChan),hmeuReco_y[(nChan+1)*iTime + iChan]);
117  }
118  }
119  }
120  }
121 
122  // Check is open
123  bool isLoaded=true;
124  int iFalse=0;
125  for( int iTime=0; iTime<=timeDiv; ++iTime ){
126  if( runsPerDiv>0 && iTime==timeDiv ) break; // include iTime==0 when that's the only one, otherwise break when iTime==timeDiv
127  if(!hmeuReco_x[(nChan+1)*iTime]){ isLoaded=false; iFalse=0;}
128  else if(!hmeuTrue_x[(nChan+1)*iTime]){ isLoaded=false; iFalse=1;}
129  else if(!hmeuReco_y[(nChan+1)*iTime]){ isLoaded=false; iFalse=2;}
130  else if(!hmeuTrue_y[(nChan+1)*iTime]){ isLoaded=false; iFalse=3;}
131  else if( byDB ){
132  for( unsigned int iChan=1; iChan<=nChan; ++iChan ){
133  if(!hmeuReco_x[(nChan+1)*iTime + iChan]){ isLoaded=false; iFalse=4;}
134  else if(!hmeuTrue_x[(nChan+1)*iTime + iChan]){ isLoaded=false; iFalse=5;}
135  else if(!hmeuReco_y[(nChan+1)*iTime + iChan]){ isLoaded=false; iFalse=6;}
136  else if(!hmeuTrue_y[(nChan+1)*iTime + iChan]){ isLoaded=false; iFalse=7;}
137  }
138  }
139  if(!isLoaded){
140  std::cout << "\nERROR!!!!! \n\t Histogram " << iFalse << " didn't load for timeDiv " << iTime << std::endl;
141  return;
142  }
143  }
144 
145  std::cout << "Writing csv ";
146 
147  for( int iTime=0; iTime<=timeDiv; ++iTime ){
148  if( runsPerDiv>0 && iTime==timeDiv ) break; // include iTime==0 when that's the only one, otherwise break when iTime==timeDiv
149  if(calibMethod.compare("nom")==0){
150  // Set the range nominally used for calibration -- slightly truncated
151  hmeuReco_x[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,100);
152  hmeuTrue_x[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
153  hmeuReco_y[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,100);
154  hmeuTrue_y[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
155  meu_x[(nChan+1)*iTime] = hmeuReco_x[(nChan+1)*iTime]->GetMean(1);
156  meu_y[(nChan+1)*iTime] = hmeuReco_y[(nChan+1)*iTime]->GetMean(1);
157  mev_x[(nChan+1)*iTime] = hmeuTrue_x[(nChan+1)*iTime]->GetMean(1);
158  mev_y[(nChan+1)*iTime] = hmeuTrue_y[(nChan+1)*iTime]->GetMean(1);
159  meuRMS_x[(nChan+1)*iTime] = hmeuReco_x[(nChan+1)*iTime]->GetRMS(1);
160  meuRMS_y[(nChan+1)*iTime] = hmeuReco_y[(nChan+1)*iTime]->GetRMS(1);
161  mevRMS_x[(nChan+1)*iTime] = hmeuTrue_x[(nChan+1)*iTime]->GetRMS(1);
162  mevRMS_y[(nChan+1)*iTime] = hmeuTrue_y[(nChan+1)*iTime]->GetRMS(1);
163  meuErr_x[(nChan+1)*iTime] = meuRMS_x[(nChan+1)*iTime]/std::sqrt(hmeuReco_x[(nChan+1)*iTime]->GetEntries());
164  meuErr_y[(nChan+1)*iTime] = meuRMS_y[(nChan+1)*iTime]/std::sqrt(hmeuReco_y[(nChan+1)*iTime]->GetEntries());
165  mevErr_x[(nChan+1)*iTime] = mevRMS_x[(nChan+1)*iTime]/std::sqrt(hmeuTrue_x[(nChan+1)*iTime]->GetEntries());
166  mevErr_y[(nChan+1)*iTime] = mevRMS_y[(nChan+1)*iTime]/std::sqrt(hmeuTrue_y[(nChan+1)*iTime]->GetEntries());
167  meu[(nChan+1)*iTime] = (meu_x[(nChan+1)*iTime]+meu_y[(nChan+1)*iTime])/2;
168  mev[(nChan+1)*iTime] = (mev_x[(nChan+1)*iTime]+mev_y[(nChan+1)*iTime])/2;
169  //dx = 1/sqrt(1/dx1^2 +1/dx^2)
170  meuErr[(nChan+1)*iTime] = 1 / std::sqrt(std::pow(meuErr_x[(nChan+1)*iTime],-2) + std::pow(meuErr_y[(nChan+1)*iTime],-2));
171  mevErr[(nChan+1)*iTime] = 1 / std::sqrt(std::pow(mevErr_x[(nChan+1)*iTime],-2) + std::pow(mevErr_y[(nChan+1)*iTime],-2));
172  // If verbosity is turned on then print out the values for the "standard" run for this time set
173  if( verbosity ){
174  std::cout << "" << std::endl;
175  std::cout << "X View: meu=" << meu_x[(nChan+1)*iTime] << ", nhits=" << hmeuReco_x[(nChan+1)*iTime]->GetEntries() << std::endl;
176  std::cout << "Y View: meu=" << meu_y[(nChan+1)*iTime] << ", nhits=" << hmeuReco_y[(nChan+1)*iTime]->GetEntries() << std::endl;
177  }
178  if( byDB ){
179  // FOR THE MC VALUES, USE THE VERSION --NOT-- BROKEN INTO DB
180  for( unsigned int iChan=1; iChan<=nChan; ++iChan ){
181  hmeuReco_x[(nChan+1)*iTime + iChan]->GetXaxis()->SetRangeUser(0,100);
182  hmeuTrue_x[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
183  hmeuReco_y[(nChan+1)*iTime + iChan]->GetXaxis()->SetRangeUser(0,100);
184  hmeuTrue_y[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
185  meu_x[(nChan+1)*iTime + iChan] = hmeuReco_x[(nChan+1)*iTime + iChan]->GetMean(1);
186  meu_y[(nChan+1)*iTime + iChan] = hmeuReco_y[(nChan+1)*iTime + iChan]->GetMean(1);
187  mev_x[(nChan+1)*iTime + iChan] = hmeuTrue_x[(nChan+1)*iTime]->GetMean(1);
188  mev_y[(nChan+1)*iTime + iChan] = hmeuTrue_y[(nChan+1)*iTime]->GetMean(1);
189  meuRMS_x[(nChan+1)*iTime + iChan] = hmeuReco_x[(nChan+1)*iTime + iChan]->GetRMS(1);
190  meuRMS_y[(nChan+1)*iTime + iChan] = hmeuReco_y[(nChan+1)*iTime + iChan]->GetRMS(1);
191  mevRMS_x[(nChan+1)*iTime + iChan] = hmeuTrue_x[(nChan+1)*iTime]->GetRMS(1);
192  mevRMS_y[(nChan+1)*iTime + iChan] = hmeuTrue_y[(nChan+1)*iTime]->GetRMS(1);
193  meuErr_x[(nChan+1)*iTime + iChan] = meuRMS_x[(nChan+1)*iTime + iChan]/std::sqrt(hmeuReco_x[(nChan+1)*iTime + iChan]->GetEntries());
194  meuErr_y[(nChan+1)*iTime + iChan] = meuRMS_y[(nChan+1)*iTime + iChan]/std::sqrt(hmeuReco_y[(nChan+1)*iTime + iChan]->GetEntries());
195  mevErr_x[(nChan+1)*iTime + iChan] = mevRMS_x[(nChan+1)*iTime + iChan]/std::sqrt(hmeuTrue_x[(nChan+1)*iTime + iChan]->GetEntries());
196  mevErr_y[(nChan+1)*iTime + iChan] = mevRMS_y[(nChan+1)*iTime + iChan]/std::sqrt(hmeuTrue_y[(nChan+1)*iTime + iChan]->GetEntries());
197  meu[(nChan+1)*iTime + iChan] = (meu_x[(nChan+1)*iTime + iChan]+meu_y[(nChan+1)*iTime + iChan])/2;
198  mev[(nChan+1)*iTime + iChan] = (mev_x[(nChan+1)*iTime + iChan]+mev_y[(nChan+1)*iTime + iChan])/2;
199  //dx = 1/sqrt(1/dx1^2 +1/dx^2)
200  meuErr[(nChan+1)*iTime + iChan] = 1 / std::sqrt(std::pow(meuErr_x[(nChan+1)*iTime + iChan],-2) + std::pow(meuErr_y[(nChan+1)*iTime + iChan],-2));
201  mevErr[(nChan+1)*iTime + iChan] = 1 / std::sqrt(std::pow(mevErr_x[(nChan+1)*iTime + iChan],-2) + std::pow(mevErr_y[(nChan+1)*iTime + iChan],-2));
202  }
203  } // end if diblock
204  } // end if nominal
205  else if(calibMethod.compare("mpv")==0){
206  // Still using RMS for error (as in nom), but use the center of
207  // max bin for Most Probable Value (MPV)
208  hmeuReco_x[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,100);
209  hmeuTrue_x[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
210  hmeuReco_y[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,100);
211  hmeuTrue_y[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
212  meu_x[(nChan+1)*iTime] = hmeuReco_x[(nChan+1)*iTime]->GetXaxis()->GetBinCenter( hmeuReco_x[(nChan+1)*iTime]->GetMaximumBin() );
213  meu_y[(nChan+1)*iTime] = hmeuReco_y[(nChan+1)*iTime]->GetXaxis()->GetBinCenter( hmeuReco_y[(nChan+1)*iTime]->GetMaximumBin() );
214  mev_x[(nChan+1)*iTime] = hmeuTrue_x[(nChan+1)*iTime]->GetXaxis()->GetBinCenter( hmeuTrue_x[(nChan+1)*iTime]->GetMaximumBin() );
215  mev_y[(nChan+1)*iTime] = hmeuTrue_y[(nChan+1)*iTime]->GetXaxis()->GetBinCenter( hmeuTrue_y[(nChan+1)*iTime]->GetMaximumBin() );
216  meuRMS_x[(nChan+1)*iTime] = hmeuReco_x[(nChan+1)*iTime]->GetRMS(1);
217  meuRMS_y[(nChan+1)*iTime] = hmeuReco_y[(nChan+1)*iTime]->GetRMS(1);
218  mevRMS_x[(nChan+1)*iTime] = hmeuTrue_x[(nChan+1)*iTime]->GetRMS(1);
219  mevRMS_y[(nChan+1)*iTime] = hmeuTrue_y[(nChan+1)*iTime]->GetRMS(1);
220  meuErr_x[(nChan+1)*iTime] = meuRMS_x[(nChan+1)*iTime]/std::sqrt(hmeuReco_x[(nChan+1)*iTime]->GetEntries());
221  meuErr_y[(nChan+1)*iTime] = meuRMS_y[(nChan+1)*iTime]/std::sqrt(hmeuReco_y[(nChan+1)*iTime]->GetEntries());
222  mevErr_x[(nChan+1)*iTime] = mevRMS_x[(nChan+1)*iTime]/std::sqrt(hmeuTrue_x[(nChan+1)*iTime]->GetEntries());
223  mevErr_y[(nChan+1)*iTime] = mevRMS_y[(nChan+1)*iTime]/std::sqrt(hmeuTrue_y[(nChan+1)*iTime]->GetEntries());
224  meu[(nChan+1)*iTime] = (meu_x[(nChan+1)*iTime]+meu_y[(nChan+1)*iTime])/2;
225  mev[(nChan+1)*iTime] = (mev_x[(nChan+1)*iTime]+mev_y[(nChan+1)*iTime])/2;
226  //dx = 1/sqrt(1/dx1^2 +1/dx^2)
227  meuErr[(nChan+1)*iTime] = 1 / std::sqrt(std::pow(meuErr_x[(nChan+1)*iTime],-2) + std::pow(meuErr_y[(nChan+1)*iTime],-2));
228  mevErr[(nChan+1)*iTime] = 1 / std::sqrt(std::pow(mevErr_x[(nChan+1)*iTime],-2) + std::pow(mevErr_y[(nChan+1)*iTime],-2));
229  if( verbosity ){
230  std::cout << "" << std::endl;
231  std::cout << "X View: meu=" << meu_x[(nChan+1)*iTime] << ", nhits=" << hmeuReco_x[(nChan+1)*iTime]->GetEntries() << std::endl;
232  std::cout << "Y View: meu=" << meu_y[(nChan+1)*iTime] << ", nhits=" << hmeuReco_y[(nChan+1)*iTime]->GetEntries() << std::endl;
233  }
234  if( byDB ){
235  // FOR THE MC VALUES, USE THE VERSION --NOT-- BROKEN INTO DB
236  for( unsigned int iChan=1; iChan<=nChan; ++iChan ){
237  hmeuReco_x[(nChan+1)*iTime + iChan]->GetXaxis()->SetRangeUser(0,100);
238  hmeuTrue_x[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
239  hmeuReco_y[(nChan+1)*iTime + iChan]->GetXaxis()->SetRangeUser(0,100);
240  hmeuTrue_y[(nChan+1)*iTime]->GetXaxis()->SetRangeUser(0,6);
241  meu_x[(nChan+1)*iTime + iChan] = hmeuReco_x[(nChan+1)*iTime + iChan]->GetXaxis()->GetBinCenter( hmeuReco_x[(nChan+1)*iTime + iChan]->GetMaximumBin() );
242  meu_y[(nChan+1)*iTime + iChan] = hmeuReco_y[(nChan+1)*iTime + iChan]->GetXaxis()->GetBinCenter( hmeuReco_y[(nChan+1)*iTime + iChan]->GetMaximumBin() );
243  mev_x[(nChan+1)*iTime + iChan] = hmeuTrue_x[(nChan+1)*iTime]->GetXaxis()->GetBinCenter( hmeuTrue_x[(nChan+1)*iTime + iChan]->GetMaximumBin() );
244  mev_y[(nChan+1)*iTime + iChan] = hmeuTrue_y[(nChan+1)*iTime]->GetXaxis()->GetBinCenter( hmeuTrue_y[(nChan+1)*iTime + iChan]->GetMaximumBin() );
245  meuRMS_x[(nChan+1)*iTime + iChan] = hmeuReco_x[(nChan+1)*iTime + iChan]->GetRMS(1);
246  meuRMS_y[(nChan+1)*iTime + iChan] = hmeuReco_y[(nChan+1)*iTime + iChan]->GetRMS(1);
247  mevRMS_x[(nChan+1)*iTime + iChan] = hmeuTrue_x[(nChan+1)*iTime]->GetRMS(1);
248  mevRMS_y[(nChan+1)*iTime + iChan] = hmeuTrue_y[(nChan+1)*iTime]->GetRMS(1);
249  meuErr_x[(nChan+1)*iTime + iChan] = meuRMS_x[(nChan+1)*iTime + iChan]/std::sqrt(hmeuReco_x[(nChan+1)*iTime + iChan]->GetEntries());
250  meuErr_y[(nChan+1)*iTime + iChan] = meuRMS_y[(nChan+1)*iTime + iChan]/std::sqrt(hmeuReco_y[(nChan+1)*iTime + iChan]->GetEntries());
251  mevErr_x[(nChan+1)*iTime + iChan] = mevRMS_x[(nChan+1)*iTime + iChan]/std::sqrt(hmeuTrue_x[(nChan+1)*iTime + iChan]->GetEntries());
252  mevErr_y[(nChan+1)*iTime + iChan] = mevRMS_y[(nChan+1)*iTime + iChan]/std::sqrt(hmeuTrue_y[(nChan+1)*iTime + iChan]->GetEntries());
253  meu[(nChan+1)*iTime + iChan] = (meu_x[(nChan+1)*iTime + iChan]+meu_y[(nChan+1)*iTime + iChan])/2;
254  mev[(nChan+1)*iTime + iChan] = (mev_x[(nChan+1)*iTime + iChan]+mev_y[(nChan+1)*iTime + iChan])/2;
255  //dx = 1/sqrt(1/dx1^2 +1/dx^2)
256  meuErr[(nChan+1)*iTime + iChan] = 1 / std::sqrt(std::pow(meuErr_x[(nChan+1)*iTime + iChan],-2) + std::pow(meuErr_y[(nChan+1)*iTime + iChan],-2));
257  mevErr[(nChan+1)*iTime + iChan] = 1 / std::sqrt(std::pow(mevErr_x[(nChan+1)*iTime + iChan],-2) + std::pow(mevErr_y[(nChan+1)*iTime + iChan],-2));
258  }
259  } // end if diblock
260  } // end if mpv
261  // NOTE: YOUR METHOD HERE!
262  // if( calibMethod.find("")!=std::string::npos ){
263  //
264  // }
265  else{
266  std::cout << "ERROR!!!!"
267  << "\n\t --calib not set to a viable option... You passed it " << calibMethod
268  << "\n\t Returning..." << std::endl;
269  return;
270  }
271 
272 
273  // Write the results to csv files
274  std::string csvName;
275  if(timeDiv==0 && metadataPlot) csvName = getOutputName(detector,period,"",outputName,true,user_dir)+"_"+calibMethod+"_r"+std::to_string(loRun)+"-r"+std::to_string(hiRun)+".csv";
276  else if(timeDiv==0 && !metadataPlot) csvName = getOutputName(detector,period,"",outputName,true,user_dir)+"_"+calibMethod+".csv";
277  else csvName = getOutputName(detector,period,"",outputName,true,user_dir)+
278  "_"+calibMethod+
279  "_r"+std::to_string(loRun+(iTime*runsPerDiv))+
280  "-r"+std::to_string(loRun+((iTime+1)*runsPerDiv)-1>hiRun ? hiRun : loRun+((iTime+1)*runsPerDiv)-1)+".csv";
281  std::cout << csvName.c_str() << " ... " << std::endl;
282  std::ofstream csv( csvName.c_str() );
283  if(!csv.is_open()){
284  std::cout << "File not found" << std::endl;
285  return;
286  }
287  csv << "#meu,meu_err,gev,channel,tv" << std::endl;
288  // meu, meu err, gev, channel, tv
289  for( unsigned int iChan=1; iChan<=nChan; ++iChan ){
290  csv << Form("%f,%f,%f,%i,%i", (!byDB ? meu[(nChan+1)*iTime] : meu[(nChan+1)*iTime + iChan]),
291  (!byDB ? meuErr[(nChan+1)*iTime] : meuErr[(nChan+1)*iTime + iChan]),
292  0.001*(!byDB ? mev[(nChan+1)*iTime] : mev[(nChan+1)*iTime + iChan]),
293  iChan,
294  tv) << std::endl;
295  }
296  csv.close();
297  } // end time div
298 
299  std::cout << "Done!" << std::endl;
300 
301  return;
302 }
303 
304 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////
305 // Make the file with plots used to perform the calibration
306 //////////////////////////////////////////////////////////////////////////////////////////////////////////////////
307 void RunCalibration(TFile *f, std::string det, bool isMC, std::string inFileName, std::string outFileName, long maxEntries, bool byDB, int runsPerDiv, bool checkNoSwap){
308  if( !f ){
309  std::cout << "ERROR!"
310  << "\n\t You must supply a valid file to read in hit info..."
311  << "\n\t Your file isn't open!"
312  << std::endl;
313  return;
314  }
315 
316  // TODO: this should work, but is clunky. Maybe better way?
317  bool ProtondEdx = false;
318  if( inFileName.find("/nova/ana/users/sgermani")!=std::string::npos ){
319  if( outFileName.find("proton")!=std::string::npos ) ProtondEdx = true;
320  }
321 
322  // NOTE: this is commented out because unused...
323  //int SelectPDG = -1; //-1 for all
324 
325  if(isMC) std::cout << "\n MC file, plot true values as well \n" << std::endl;
326  else std::cout << "\n Data file, no true values to plot \n" << std::endl;
327 
328  unsigned int nChan = 14;
329  if( det.find("fd")!=std::string::npos ) nChan = 14;
330  else if( det.find("nd")!=std::string::npos ) nChan = 4;
331 
332  TFile *outFile = new TFile(outFileName.c_str(),"recreate");
333 
334  TTree *t;
335  std::vector< std::string > forest; forest.resize(4);
336  forest[0] = "muondedxana/fTree";
337  forest[1] = "remiddedx/fDedxSample";
338  forest[2] = "assess/tree";
339  forest[3] = "remiddedxfd/fDedxSample";
340  if(ProtondEdx)
341  forest[1] = "remiddedx/fProtonDedxSample";
342  unsigned int tree_i=999;
343  for(unsigned int i=0; i<forest.size(); i++){
344  f->GetObject(forest[i].c_str(),t);
345  tree_i = i;
346  if(t) break;
347  }
348  if(!t) std::cout << "could not find tree" << std::endl;
349 
350  float PECorr, PE, x, w, true_path, path, trueE, recoGeV;
351  float xDir, yDir, zDir, dcosz, dcosy;
352  float fbBin_fl, cell_fl, plane_fl;
353  int cell, plane, diblock, fbBin, view, run;
354  int nXYHits;
355  float xpos, ypos, adc;
356  int pdg, cells_per_plane;
357  uint32_t EventTimeHigh;
358 
359  if( tree_i==2 ){
360  t->SetBranchAddress( "cell", &cell );
361  t->SetBranchAddress( "plane", &plane );
362  t->SetBranchAddress( "W", &w );
363  t->SetBranchAddress( "view", &view );
364  //t->SetBranchAddress( "path", &path );
365  t->SetBranchAddress( "pecorr", &PECorr );
366  t->SetBranchAddress( "pe", &PE );
367  t->SetBranchAddress( "trueE", &trueE ); // gev
368  t->SetBranchAddress( "fibBrightnessBin",&fbBin);
369  x=0;
370  true_path = 0;
371  path = 0;
372 
373  } else if (tree_i==1 || tree_i==3){
374 
375  if(ProtondEdx){
376  t->SetBranchAddress( "Protoncell", &cell );
377  t->SetBranchAddress( "Protonplane", &plane );
378  t->SetBranchAddress( "ProtondistFromEnd", &x );
379  t->SetBranchAddress( "ProtonxHit", &xpos );
380  t->SetBranchAddress( "ProtonyHit", &ypos );
381  t->SetBranchAddress( "Protonview", &view );
382  t->SetBranchAddress( "ProtonxDir", &xDir );
383  t->SetBranchAddress( "ProtonyDir", &yDir );
384  t->SetBranchAddress( "ProtonzDir", &zDir );
385  t->SetBranchAddress( "Protondx", &path );
386  t->SetBranchAddress( "Protonmcdx", &true_path );
387  t->SetBranchAddress( "Protonpecorr", &PECorr );
388  t->SetBranchAddress( "Protonpe", &PE );
389  t->SetBranchAddress( "Protonmcde", &trueE ); // gev
390  t->SetBranchAddress( "Protonmcdx", &true_path ); // cm
391  t->SetBranchAddress( "ProtonfiberBrightness",&fbBin);
392  t->SetBranchAddress( "Protonpdg", &pdg );
393  t->SetBranchAddress( "Protonde", &recoGeV );
394  t->SetBranchAddress( "Protonadc", &adc );
395  t->SetBranchAddress( "ProtonncellsPerPlane", &cells_per_plane );
396  } else {
397  t->SetBranchAddress( "cell", &cell );
398  t->SetBranchAddress( "plane", &plane );
399  // t->SetBranchAddress( "diblock", &diblock );
400  // t->SetBranchAddress( "run", &run );
401  t->SetBranchAddress( "distFromEnd", &x );
402  t->SetBranchAddress( "xHit", &xpos );
403  t->SetBranchAddress( "yHit", &ypos );
404  t->SetBranchAddress( "view", &view );
405  t->SetBranchAddress( "xDir", &xDir );
406  t->SetBranchAddress( "yDir", &yDir );
407  t->SetBranchAddress( "zDir", &zDir );
408  t->SetBranchAddress( "dx", &path );
409  t->SetBranchAddress( "mcdx", &true_path );
410  t->SetBranchAddress( "pecorr", &PECorr );
411  t->SetBranchAddress( "pe", &PE );
412  //t->SetBranchAddress( "trueE", &trueE ); // this was true E of MCParticle (GeV)
413  t->SetBranchAddress( "mcde", &trueE ); // gev
414  t->SetBranchAddress( "mcdx", &true_path ); // cm
415  t->SetBranchAddress( "fiberBrightness",&fbBin);
416  t->SetBranchAddress( "pdg", &pdg );
417  t->SetBranchAddress( "de", &recoGeV );
418  t->SetBranchAddress( "adc", &adc );
419  t->SetBranchAddress( "ncellsPerPlane", &cells_per_plane );
420  }
421 
422  } else {
423  t->SetBranchAddress( "nXYHits", &nXYHits); // number of hits (I think per track). Should weight by 1/nXYHits to end up getting 1 per track then...
424  t->SetBranchAddress( "x", &x );
425  t->SetBranchAddress( "pdg", &pdg );
426  t->SetBranchAddress( "w", &w );
427  t->SetBranchAddress( "view", &view );
428  t->SetBranchAddress( "path", &path );
429  t->SetBranchAddress( "truepath", &true_path );
430  t->SetBranchAddress( "cell", &cell_fl );
431  t->SetBranchAddress( "plane", &plane_fl );
432  t->SetBranchAddress( "diblock", &diblock ); // the other stuff is _fl but I'm pretty sure int is the right thing for this
433  t->SetBranchAddress( "run", &run ); // the other stuff is _fl but I'm pretty sure int is the right thing for this
434  t->SetBranchAddress( "EventTimeHigh", &EventTimeHigh ); // use the timestamp itself w/o needing to convert the evt_ times
435  t->SetBranchAddress( "PECorr", &PECorr );
436  t->SetBranchAddress( "PE", &PE );
437  t->SetBranchAddress( "trueE", &trueE ); // mev
438  t->SetBranchAddress( "fibBrightnessBin",&fbBin_fl);
439  t->SetBranchAddress( "RecoGeV", &recoGeV );
440  t->SetBranchAddress( "dCosZ", &dcosz );
441  t->SetBranchAddress( "dCosY", &dcosy );
442  t->SetBranchAddress( "diblock", &diblock );
443  }
444 
445 
446  ///////////////////////////////////////////////////////////////////
447  ///////////////////////////////////////////////////////////////////
448  ////////////////// //////////////////
449  ////////////////// Set up plot ranges //////////////////
450  ////////////////// //////////////////
451  ///////////////////////////////////////////////////////////////////
452  ///////////////////////////////////////////////////////////////////
453  unsigned int nFBbins = 9;
454  int nCells(0), PlanesPerView(0), tprof_bins(30), dxbins(30);
455  float wrange;
456  float dxmin(3), dxmax(10);
457  float pe_x_min, pe_x_max, pe_y_min, pe_y_max;
458  float track_end_max(500); // for plotting regardless of track window
459  // NOTE: this is commented out because unused...
460  //float w_x_cut_min, w_x_cut_max, w_y_cut_min, w_y_cut_max;
461  float w_y_cut_min, w_y_cut_max;
462  float flatW_x_min, flatW_x_max, flatW_y_min, flatW_y_max, flatW_x_selRgn, flatW_y_selRgn;
463  if(det=="nd"){
464 
465  nCells = 32*3;
466  PlanesPerView = 96;
467  wrange = 200;
468 
469  pe_x_min = 0; pe_x_max = 450;
470  pe_y_min = 0; pe_y_max = 450;
471  // NOTE: this is commented out because unused...
472  //w_x_cut_min = -100; w_x_cut_max = 150;
473  w_y_cut_min = -150; w_y_cut_max = 150;
474 
475  flatW_x_min = -100; flatW_x_max = 100;
476  flatW_y_min = -100; flatW_y_max = 100;
477 
478  flatW_x_selRgn = (1./3.)*(flatW_x_max-flatW_x_min);
479  flatW_y_selRgn = (1./3.)*(flatW_y_max-flatW_y_min);
480 
481  //if(isMC) wbins = 15;
482 
483  if(ProtondEdx){
484  tprof_bins=15;
485  dxbins=30;
486  pe_y_max = 1300;
487  pe_x_max = 1300;
488  track_end_max = 500;
489  }
490 
491 
492  } else {
493 
494  nCells = 384;
495  PlanesPerView = 448;
496  wrange = 800;
497 
498  pe_x_min = 0; pe_x_max = 450;
499  pe_y_min = 0; pe_y_max = 450;
500  // NOTE: this is commented out because unused...
501  //w_x_cut_min = -600; w_x_cut_max = 600;
502  w_y_cut_min = -600; w_y_cut_max = 600;
503 
504  flatW_x_min = 200; flatW_x_max = 600;
505  flatW_y_min = 200; flatW_y_max = 600;
506 
507  flatW_x_selRgn = (1./3.)*(flatW_x_max-flatW_x_min);
508  flatW_y_selRgn = (1./3.)*(flatW_y_max-flatW_y_min);
509  }
510 
511 
512  std::cout << "define plots" << std::endl;
513 
514  //////////////////////////////////////////////////////
515  ///////////// METADATA PLOT /////////////
516  // -------------------- NEW!!!! --------------------
517  if( (byDB || runsPerDiv!=0) && tree_i!=0 ){
518  std::cout << "ERROR!!!!!!!"
519  << "\n\t Not set up to byDB or run divisions with this type of calib tree... Returning."
520  << std::endl;
521  return;
522  }
523 
524  unsigned int nDivisions = 0;
525  int loRun(0), hiRun(0);
526  // Run through the events once to determine the min and max run numbers
527  // TODO: better to do this or t->GetMaximum("run") ???
528  unsigned int ncheck = t->GetEntries();
529  unsigned int nskim = ncheck;
530  if( maxEntries>0 && ncheck>maxEntries ) nskim = maxEntries;
531  for( unsigned int iEntry=0; iEntry<nskim; ++iEntry ){
532  t->GetEntry(iEntry);
533  if( iEntry==0 ){
534  loRun = run;
535  hiRun = run;
536  }
537  else{
538  if( run<loRun ) loRun = run;
539  if( run>hiRun ) hiRun = run;
540  }
541  } // end find lo and hi run numbers
542 
543  if( runsPerDiv!=0 ){
544  nDivisions = ceil( float(hiRun-loRun+1)/float(runsPerDiv) );
545  } // end set nDivisions
546 
547  if( runsPerDiv!=0 && nDivisions==0 ){
548  std::cout << "ERROR!!!!!!!"
549  << "\n\t Something is wrong...runsPerDiv!=0 implies you want run divisions, but nDivisions=0..."
550  << "\n\t Returning."
551  << std::endl;
552  return;
553  }
554 
555  std::cout << "\t" << "loRun = " << loRun << ", and hiRun = " << hiRun << std::endl;
556  std::cout << "\t" << "nDivisions = " << nDivisions << " and byDB = " << byDB << std::endl;
557 
558  TH1D *metadataPlot = new TH1D("metadataPlot",";byDB,runsPerDiv,loRun,hiRun,nDivisions;Value",5,0,5);
559  metadataPlot->SetBinContent(1,(byDB ? 1. : 0.));
560  metadataPlot->SetBinContent(2,float(runsPerDiv)); // 0 means no run breakdown
561  metadataPlot->SetBinContent(3,float(loRun)); // low run in the file
562  metadataPlot->SetBinContent(4,float(hiRun)); // high run in the file
563  metadataPlot->SetBinContent(5,float(nDivisions)); // 0 means no run breakdown
564  //////////////////////////////////////////////////////
565 
566  std::cout << " metadata plot made." << std::endl;
567 
568 
569  ///////////////////////////////////////////////////////////////////
570  ///////////////////////////////////////////////////////////////////
571  ////////////////// //////////////////
572  ////////////////// Declare plots with strict //////////////////
573  ////////////////// naming convention: //////////////////
574  /////// ///////
575  /////// "plottype_fillweight_yaxis(-xaxis)_specifiers" ///////
576  /////// ///////
577  /////// root object name is same as ///////
578  /////// identifier, except with "-" becoming "_" ///////
579  /////// ///////
580  ///////////////////////////////////////////////////////////////////
581  ///////////////////////////////////////////////////////////////////
582 
583  // updated name convention: replace just "X view" or "Y view" in title with
584  // "<Y axis> vs. <X axis> (X/Y view)"
585  // in all cases where possible. In case of TH2D, where possible, do
586  // "<Y axis> vs. <X axis> (Z axis, X/Y view)"
587 
588  TH2D *th2d_nhits_dE_pe_xview = new TH2D("th2d_nhits_dE-pe_xview",
589  "True Deposition vs. PE (n Hits, X view);PE;True Deposition (MeV)",
590  tprof_bins,0,300,
591  tprof_bins,0,0.09*1000);
592  TH2D *th2d_nhits_dE_pe_yview = new TH2D("th2d_nhits_dE-pe_yview",
593  "True Deposition vs. PE (n Hits, Y view);PE;True Deposition (MeV)",
594  tprof_bins,0,300,
595  tprof_bins,0,0.09*1000);
596 
597 
598  /* # hits
599 
600  th2d_nhits_cell-w_xview
601  th2d_nhits_w-cell_yview
602 
603  */
604 
605 
606  TH1D *th1d_nhits_resp_flatW = new TH1D("th1d_nhits_resp_flatW","Number of hits vs. PECorr/cm;PECorr/cm;# Hits",
607  1000,0,500);
608  TH1D *th1d_nhits_truededx_flatW = new TH1D("th1d_nhits_truededx_flatW","Number of hits vs. True dE/dx;MeV/cm;# Hits",
609  1000,0,25);
610  TH1D *th1d_nhits_resp_flatW_xview = new TH1D("th1d_nhits_resp_flatW_xview","Number of hits vs. PECorr/cm (X view);PECorr/cm;# Hits",
611  1000,0,500);
612  TH1D *th1d_nhits_resp_flatW_yview = new TH1D("th1d_nhits_resp_flatW_yview","Number of hits vs. PECorr/cm (Y view);PECorr/cm;# Hits",
613  1000,0,500);
614  TH1D *th1d_nhits_truededx_flatW_xview = new TH1D("th1d_nhits_truededx_flatW_xview",
615  "Number of hits vs. True dE/dx (X view);True MeV/cm;# Hits",
616  1000,0,25);
617  TH1D *th1d_nhits_truededx_flatW_yview = new TH1D("th1d_nhits_truededx_flatW_yview",
618  "Number of hits vs. True dE/dx (Y view);True MeV/cm;# Hits",
619  1000,0,25);
620  TH1D *th1d_nhits_recodedx_flatW_xview = new TH1D("th1d_nhits_recodedx_flatW_xview",
621  "Number of hits vs. Reco dE/dx (X view);Reco MeV/cm;# Hits",
622  tprof_bins,0,6);
623  TH1D *th1d_nhits_recodedx_flatW_yview = new TH1D("th1d_nhits_recodedx_flatW_yview",
624  "Number of hits vs. Reco dE/dx (Y view);Reco MeV/cm;# Hits",
625  tprof_bins,0,6);
626 
627  // NEW!!!!!
628  // Histogram for nXYHits: Hits/track (when we fill, weight by 1/(Hits/track) to get 1 for each track)
629  TH1D *th1d_nxyhits_xview = new TH1D("th1d_nxhits_xview","Histogram of number of tricell hits (X view);nXYHits;# Occurrences",
630  1000,0,1000);
631  TH1D *th1d_nxyhits_yview = new TH1D("th1d_nxhits_yview","Histogram of number of tricell hits (Y view);nXYHits;# Occurrences",
632  1000,0,1000);
633  TH1D *th1d_nxyhits_flatW_xview = new TH1D("th1d_nxhits_flatW_xview","Histogram of number of tricell hits (X view);nXYHits;# Occurrences",
634  1000,0,1000);
635  TH1D *th1d_nxyhits_flatW_yview = new TH1D("th1d_nxhits_flatW_yview","Histogram of number of tricell hits (Y view);nXYHits;# Occurrences",
636  1000,0,1000);
637 
638  //////////////////////////////////////////////////////
639  // -------------------- NEW!!!! --------------------
640  // Plots broken down into time period and diblocks (similar naming convention to fbsplit for use in the Validation Webpage infrastructure)
641  std::vector<TH1D*> th1d_nhits_resp_flatW_diblock; // split by diblock
642  th1d_nhits_resp_flatW_diblock.resize(nChan);
643  std::vector<TH1D*> th1d_nhits_resp_flatW_xview_diblock; // split by diblock
644  th1d_nhits_resp_flatW_xview_diblock.resize(nChan);
645  std::vector<TH1D*> th1d_nhits_resp_flatW_yview_diblock; // split by diblock
646  th1d_nhits_resp_flatW_yview_diblock.resize(nChan);
647 
648  std::vector<TH1D*> th1d_nhits_truededx_flatW_diblock;
649  th1d_nhits_truededx_flatW_diblock.resize(nChan);
650  std::vector<TH1D*> th1d_nhits_truededx_flatW_xview_diblock;
651  th1d_nhits_truededx_flatW_xview_diblock.resize(nChan);
652  std::vector<TH1D*> th1d_nhits_truededx_flatW_yview_diblock;
653  th1d_nhits_truededx_flatW_yview_diblock.resize(nChan);
654 
655  std::vector<TH1D*> th1d_nhits_resp_flatW_range; // split by run division
656  if( nDivisions>0 ) th1d_nhits_resp_flatW_range.resize(nDivisions);
657  std::vector<TH1D*> th1d_nhits_resp_flatW_xview_range; // split by run division
658  if( nDivisions>0 ) th1d_nhits_resp_flatW_xview_range.resize(nDivisions);
659  std::vector<TH1D*> th1d_nhits_resp_flatW_yview_range; // split by run division
660  if( nDivisions>0 ) th1d_nhits_resp_flatW_yview_range.resize(nDivisions);
661 
662  std::vector< std::vector<TH1D*> > th1d_nhits_resp_flatW_range_diblock; // split by range then diblock
663  if( nDivisions>0 ){
664  th1d_nhits_resp_flatW_range_diblock.resize(nDivisions);
665  for( long thisTime=0; thisTime<nDivisions; ++thisTime ){
666  th1d_nhits_resp_flatW_range_diblock[thisTime].resize(nChan);
667  }
668  }
669  std::vector< std::vector<TH1D*> > th1d_nhits_resp_flatW_xview_range_diblock; // split by range then diblock
670  if( nDivisions>0 ){
671  th1d_nhits_resp_flatW_xview_range_diblock.resize(nDivisions);
672  for( long thisTime=0; thisTime<nDivisions; ++thisTime ){
673  th1d_nhits_resp_flatW_xview_range_diblock[thisTime].resize(nChan);
674  }
675  }
676  std::vector< std::vector<TH1D*> > th1d_nhits_resp_flatW_yview_range_diblock; // split by range then diblock
677  if( nDivisions>0 ){
678  th1d_nhits_resp_flatW_yview_range_diblock.resize(nDivisions);
679  for( unsigned int thisTime=0; thisTime<nDivisions; ++thisTime ){
680  th1d_nhits_resp_flatW_yview_range_diblock[thisTime].resize(nChan);
681  }
682  }
683 
684  if( nDivisions>0 ){
685  for( unsigned int thisTime=0; thisTime<nDivisions; ++thisTime ){
686  th1d_nhits_resp_flatW_range[thisTime] = new TH1D( TString::Format("th1d_range%u_nhits_resp_flatW",thisTime),
687  "Number of hits vs. PECorr/cm;PECorr/cm;# Hits",1000,0,500);
688  th1d_nhits_resp_flatW_xview_range[thisTime] = new TH1D( TString::Format("th1d_range%u_nhits_resp_flatW_xview",thisTime),
689  "Number of hits vs. PECorr/cm (X view);PECorr/cm;# Hits",1000,0,500);
690  th1d_nhits_resp_flatW_yview_range[thisTime] = new TH1D( TString::Format("th1d_range%u_nhits_resp_flatW_yview",thisTime),
691  "Number of hits vs. PECorr/cm (Y view);PECorr/cm;# Hits",1000,0,500);
692  }
693  }
694  if( byDB ){
695  for( unsigned int iChan=0; iChan<nChan; ++iChan ){
696  th1d_nhits_resp_flatW_diblock[iChan] = new TH1D( TString::Format("th1d_db%u_nhits_resp_flatW",iChan+1),
697  "Number of hits vs. PECorr/cm;PECorr/cm;# Hits",1000,0,500);
698  th1d_nhits_resp_flatW_xview_diblock[iChan] = new TH1D( TString::Format("th1d_db%u_nhits_resp_flatW_xview",iChan+1),
699  "Number of hits vs. PECorr/cm (X view);PECorr/cm;# Hits",1000,0,500);
700  th1d_nhits_resp_flatW_yview_diblock[iChan] = new TH1D( TString::Format("th1d_db%u_nhits_resp_flatW_yview",iChan+1),
701  "Number of hits vs. PECorr/cm (Y view);PECorr/cm;# Hits",1000,0,500);
702  th1d_nhits_truededx_flatW_diblock[iChan] = new TH1D( TString::Format("th1d_db%u_nhits_truededx_flatW",iChan+1),
703  "Number of hits vs. True dE/dx;True MeV/cm;# Hits",1000,0,25);
704  th1d_nhits_truededx_flatW_xview_diblock[iChan] = new TH1D( TString::Format("th1d_db%u_nhits_truededx_flatW_xview",iChan+1),
705  "Number of hits vs. True dE/dx (X view);True MeV/cm;# Hits",1000,0,25);
706  th1d_nhits_truededx_flatW_yview_diblock[iChan] = new TH1D( TString::Format("th1d_db%u_nhits_truededx_flatW_yview",iChan+1),
707  "Number of hits vs. True dE/dx (Y view);True MeV/cm;# Hits",1000,0,25);
708  if( nDivisions>0 ){
709  for( unsigned int thisTime=0; thisTime<nDivisions; ++thisTime ){
710  th1d_nhits_resp_flatW_range_diblock[thisTime][iChan] = new TH1D( TString::Format("th1d_range%u_db%u_nhits_resp_flatW",thisTime,iChan+1),
711  "Number of hits vs. PECorr/cm;PECorr/cm;# Hits",1000,0,500);
712  th1d_nhits_resp_flatW_xview_range_diblock[thisTime][iChan] = new TH1D( TString::Format("th1d_range%u_db%u_nhits_resp_flatW_xview",thisTime,iChan+1),
713  "Number of hits vs. PECorr/cm (X view);PECorr/cm;# Hits",1000,0,500);
714  th1d_nhits_resp_flatW_yview_range_diblock[thisTime][iChan] = new TH1D( TString::Format("th1d_range%u_db%u_nhits_resp_flatW_yview",thisTime,iChan+1),
715  "Number of hits vs. PECorr/cm (Y view);PECorr/cm;# Hits",1000,0,500);
716  }
717  }
718  }
719  }
720  //////////////////////////////////////////////////////
721 
722 
723  //////////////////////////////////////////////////////
724  ///// NEW Profile plots over time
725  TProfile *tprof_resp_time_flatW_xview = new TProfile("tprof_resp_time_flatW_xview","Response [PECorr/cm] vs Time of Event;Time;Response",
726  6312,1388600000,1704100000);
727  TProfile *tprof_pecm_time_flatW_xview = new TProfile("tprof_pecm_time_flatW_xview","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
728  6312,1388600000,1704100000);
729  TProfile *tprof_pe_time_flatW_xview = new TProfile("tprof_pe_time_flatW_xview","Uncorrected PE vs Time of Event;Time;Response",
730  6312,1388600000,1704100000);
731  TProfile *tprof_resp_time_flatW_yview = new TProfile("tprof_resp_time_flatW_yview","Response [PECorr/cm] vs Time of Event;Time;Response",
732  6312,1388600000,1704100000);
733  TProfile *tprof_pecm_time_flatW_yview = new TProfile("tprof_pecm_time_flatW_yview","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
734  6312,1388600000,1704100000);
735  TProfile *tprof_pe_time_flatW_yview = new TProfile("tprof_pe_time_flatW_yview","Uncorrected PE vs Time of Event;Time;Response",
736  6312,1388600000,1704100000);
737  std::vector<TProfile*> tprof_pecm_time_flatW_xview_diblock; tprof_pecm_time_flatW_xview_diblock.resize(nChan);
738  std::vector<TProfile*> tprof_pecm_time_flatW_yview_diblock; tprof_pecm_time_flatW_yview_diblock.resize(nChan);
739  if( byDB ){
740  for( unsigned int iChan=0; iChan<tprof_pecm_time_flatW_xview_diblock.size(); ++iChan ){
741  tprof_pecm_time_flatW_xview_diblock[iChan] = new TProfile(TString::Format("tprof_pecm_time_flatW_xview_db%u",iChan),"Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
742  6312,1388600000,1704100000);
743  tprof_pecm_time_flatW_yview_diblock[iChan] = new TProfile(TString::Format("tprof_pecm_time_flatW_yview_db%u",iChan),"Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
744  6312,1388600000,1704100000);
745  }
746  }
747 
748  // Versions that cut out regions of W or channels (specifically fill with channels that have never been swapped at FD)
749  TProfile *tprof_resp_time_flatW_xview_loW = new TProfile("tprof_resp_time_flatW_xview_loW","Response [PECorr/cm] vs Time of Event;Time;Response",
750  6312,1388600000,1704100000);
751  TProfile *tprof_pecm_time_flatW_xview_loW = new TProfile("tprof_pecm_time_flatW_xview_loW","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
752  6312,1388600000,1704100000);
753  TProfile *tprof_resp_time_flatW_yview_loW = new TProfile("tprof_resp_time_flatW_yview_loW","Response [PECorr/cm] vs Time of Event;Time;Response",
754  6312,1388600000,1704100000);
755  TProfile *tprof_pecm_time_flatW_yview_loW = new TProfile("tprof_pecm_time_flatW_yview_loW","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
756  6312,1388600000,1704100000);
757  TProfile *tprof_resp_time_flatW_xview_mdW = new TProfile("tprof_resp_time_flatW_xview_mdW","Response [PECorr/cm] vs Time of Event;Time;Response",
758  6312,1388600000,1704100000);
759  TProfile *tprof_pecm_time_flatW_xview_mdW = new TProfile("tprof_pecm_time_flatW_xview_mdW","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
760  6312,1388600000,1704100000);
761  TProfile *tprof_resp_time_flatW_yview_mdW = new TProfile("tprof_resp_time_flatW_yview_mdW","Response [PECorr/cm] vs Time of Event;Time;Response",
762  6312,1388600000,1704100000);
763  TProfile *tprof_pecm_time_flatW_yview_mdW = new TProfile("tprof_pecm_time_flatW_yview_mdW","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
764  6312,1388600000,1704100000);
765  TProfile *tprof_resp_time_flatW_xview_hiW = new TProfile("tprof_resp_time_flatW_xview_hiW","Response [PECorr/cm] vs Time of Event;Time;Response",
766  6312,1388600000,1704100000);
767  TProfile *tprof_pecm_time_flatW_xview_hiW = new TProfile("tprof_pecm_time_flatW_xview_hiW","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
768  6312,1388600000,1704100000);
769  TProfile *tprof_resp_time_flatW_yview_hiW = new TProfile("tprof_resp_time_flatW_yview_hiW","Response [PECorr/cm] vs Time of Event;Time;Response",
770  6312,1388600000,1704100000);
771  TProfile *tprof_pecm_time_flatW_yview_hiW = new TProfile("tprof_pecm_time_flatW_yview_hiW","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
772  6312,1388600000,1704100000);
773 
774  TProfile *tprof_resp_time_flatW_xview_noSwap = new TProfile("tprof_resp_time_flatW_xview_noSwap","Response [PECorr/cm] vs Time of Event;Time;Response",
775  6312,1388600000,1704100000);
776  TProfile *tprof_pecm_time_flatW_xview_noSwap = new TProfile("tprof_pecm_time_flatW_xview_noSwap","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
777  6312,1388600000,1704100000);
778  TProfile *tprof_resp_time_flatW_yview_noSwap = new TProfile("tprof_resp_time_flatW_yview_noSwap","Response [PECorr/cm] vs Time of Event;Time;Response",
779  6312,1388600000,1704100000);
780  TProfile *tprof_pecm_time_flatW_yview_noSwap = new TProfile("tprof_pecm_time_flatW_yview_noSwap","Uncorrected Response [PE/cm] vs Time of Event;Time;Response",
781  6312,1388600000,1704100000);
782  //////////////////////////////////////////////////////
783 
784 
785  TH1D *th1d_nhits_truedx_xview = new TH1D("th1d_nhits_truedx_xview",
786  "Number of hits vs. True Path Length (X view);True Path Length (cm);# Hits",
787  tprof_bins,0,10);
788  TH1D *th1d_nhits_truedx_yview = new TH1D("th1d_nhits_truedx_yview",
789  "Number of hits vs. True Path Length (Y view);True Path Length (cm);# Hits",
790  tprof_bins,0,10);
791  TH1D *th1d_nhits_dx_xview = new TH1D("th1d_nhits_dx_xview",
792  "Number of hits vs. Reco Path Length (X view);Reconstructed Path Length (cm);# Hits",
793  tprof_bins,0,15);
794  TH1D *th1d_nhits_dx_yview = new TH1D("th1d_nhits_dx_yview",
795  "Number of hits vs. Reco Path Length (Y view);Reconstructed Path Length (cm);# Hits",
796  tprof_bins,0,15);
797  TH1D *th1d_nhits_dE_xview = new TH1D("th1d_nhits_dE_xview",
798  "Number of hits vs. True Energy Deposit (X view);True Energy Deposit (MeV);# Hits",
799  tprof_bins,0,20);
800  TH1D *th1d_nhits_dE_yview = new TH1D("th1d_nhits_dE_yview",
801  "Number of hits vs. True Energy Deposit (Y view);True Energy Deposit (MeV);# Hits",
802  tprof_bins,0,20);
803 
804 
805  TH1D *th1d_nhits_dcosz_xview = new TH1D("th1d_nhits_dcosz_xview",
806  "Number of Hits vs. Z Direction Cosine (X view);Z Direcion Cosine;# Hits",
807  tprof_bins,-1,1);
808  TH1D *th1d_nhits_dcosz_yview = new TH1D("th1d_nhits_dcosz_yview",
809  "Number of Hits vs. Z Direction Cosine (Y view);Z Direction Cosine;# Hits",
810  tprof_bins,-1,1);
811  TH1D *th1d_nhits_dcosy_xview = new TH1D("th1d_nhits_dcosy_xview",
812  "Number of Hits vs. Y Direction Cosine (X view);Y Direcion Cosine;# Hits",
813  tprof_bins,-1,1);
814  TH1D *th1d_nhits_dcosy_yview = new TH1D("th1d_nhits_dcosy_yview",
815  "Number of Hits vs. Y Direction Cosine (Y view);Y Direction Cosine;# Hits",
816  tprof_bins,-1,1);
817 
818 
819  TH2D *th2d_nhits_cell_w_xview = new TH2D("th2d_nhits_cell-w_xview",
820  "W vs. X Cell (n Hits, X view);X Cell;W (cm)",
821  nCells, 0, nCells,
822  tprof_bins, -wrange, wrange );
823  TH2D *th2d_nhits_w_cell_yview = new TH2D("th2d_nhits_w-cell_yview",
824  "Y Cell vs. W (n Hits, Y view);W (cm);Y Cell",
825  tprof_bins, -wrange, wrange,
826  nCells, 0, nCells );
827 
828  TH2D *th2d_nhits_dx_recodedx_xview = new TH2D("th2d_nhits_dx-recodedx_xview",
829  "Reco dE/dx vs. Reco Path Length (n Hits, X view);Reconstructed Path Length (cm);Reco dE/dx (MeV/cm)",
830  dxbins, 3, 10,
831  tprof_bins, 0, 5 );
832  TH2D *th2d_nhits_dx_recodedx_yview = new TH2D("th2d_nhits_dx-recodedx_yview",
833  "Reco dE/dx vs. Reco Path Length (n Hits, Y view);Reconstructed Path Length (cm);Reco dE/dx (MeV/cm)",
834  dxbins, 3, 10,
835  tprof_bins, 0, 5 );
836  TH2D *th2d_nhits_dx_truededx_xview = new TH2D("th2d_nhits_dx-truededx_xview",
837  "True dE/dx vs. Reco Path Length (n Hits, X view);Reconstructed Path Length (cm);True dE/dx (MeV/cm)",
838  dxbins, 3, 10,
839  tprof_bins, 0, 5 );
840  TH2D *th2d_nhits_dx_truededx_yview = new TH2D("th2d_nhits_dx-truededx_yview",
841  "True dE/dx vs. Reco Path Length (n Hits, Y view);Reconstructed Path Length (cm);True dE/dx (MeV/cm)",
842  dxbins, 3, 10,
843  tprof_bins, 0, 5 );
844 
845 
846  /* 1D profiles
847 
848  N.B. -- If there exists a fb-split vector for any particular 1d plot, the
849  combined 1d plot is kept at the end of that vector instead of here.
850 
851  tprof_resp_fbbin*
852  tprof_pecm_fbbin*
853  tprof_resp_dx_*
854  tprof_resp_cell_*
855  tprof_resp_plane_*
856  tprof_recotruedx_truedx_*
857  tprof_dxres_truedx_*
858  tprof_dxres_w_*
859 
860  */
861  TProfile *tprof_resp_fbbin = new TProfile("tprof_resp_fbbin",
862  "Response vs Fiber Brightness Bin;Fiber Brightness Bin;PECorr/cm",
863  nFBbins,0,nFBbins);
864  TProfile *tprof_pecm_fbbin = new TProfile("tprof_pecm_fbbin",
865  "Uncorrected Response vs Fiber Brightness Bin;Fiber Brightness Bin;PE/cm",
866  nFBbins,0,nFBbins);
867  TProfile *tprof_resp_fbbin_xview = new TProfile("tprof_resp_fbbin_xview",
868  "X view Response vs FB Bin;Fiber Brightness Bin;PECorr/cm",
869  nFBbins,0,nFBbins);
870  TProfile *tprof_pecm_fbbin_xview = new TProfile("tprof_pecm_fbbin_xview",
871  "X view Uncorrected Response vs FB Bin;Fiber Brightness Bin;PE/cm",
872  nFBbins,0,nFBbins);
873  TProfile *tprof_resp_fbbin_yview = new TProfile("tprof_resp_fbbin_yview",
874  "Y view Response vs FB Bin;Fiber Brightness Bin;PECorr/cm",
875  nFBbins,0,nFBbins);
876  TProfile *tprof_pecm_fbbin_yview = new TProfile("tprof_pecm_fbbin_yview",
877  "Y view Uncorrected Response vs FB Bin;Fiber Brightness Bin;PE/cm",
878  nFBbins,0,nFBbins);
879 
880  TProfile *tprof_resp_dx_xview = new TProfile("tprof_resp_dx_xview",
881  "Mean Response (X view);Path Length (cm);PECorr/cm",
882  dxbins,dxmin,dxmax);
883  TProfile *tprof_resp_dx_yview = new TProfile("tprof_resp_dx_yview",
884  "Mean Response (Y view);Path Length (cm);PECorr/cm",
885  dxbins,dxmin,dxmax);
886 
887  TProfile* tprof_pecm_cell_xview = new TProfile("tprof_pecm_cell_xview",
888  "Uncorrected Response vs. Cell Number (X view);Cell Number;PE/cm",
889  nCells, 0, nCells); // should I be super consistent when I use 'uncorrected' and such?
890  TProfile* tprof_pecm_cell_yview = new TProfile("tprof_pecm_cell_yview",
891  "Uncorrected Response vs. Cell Number (Y view);Cell Number;PE/cm",
892  nCells, 0, nCells); // should I be super consistent when I use 'uncorrected' and such?
893 
894  TProfile* tprof_resp_cell_xview = new TProfile("tprof_resp_cell_xview",
895  "Response vs Cell (X view);Cell Number;Reco Response (PECorr/cm)",
896  nCells, 0, nCells);
897  TProfile* tprof_resp_cell_yview = new TProfile("tprof_resp_cell_yview",
898  "Response vs Cell (Y view);Cell Number;Reco Response (PECorr/cm)",
899  nCells, 0, nCells);
900 
901  TProfile* tprof_resp_plane_xview = new TProfile("tprof_resp_plane_xview",
902  "Response vs Plane (X view);Plane Number;Reco Response (PECorr/cm)",
903  PlanesPerView, 0, 2*PlanesPerView);
904  TProfile* tprof_resp_plane_yview = new TProfile("tprof_resp_plane_yview",
905  "Response vs Plane (Y view);Plane Number;Reco Response (PECorr/cm)",
906  PlanesPerView, 0, 2*PlanesPerView);
907 
908  TProfile* tprof_recotruedx_truedx_xview = new TProfile("tprof_recotruedx_truedx_xview",
909  "Reco/True PCHit Path Length (X View);True Path Length (cm);Reco/True Path Length",
910  tprof_bins,1,10);
911  TProfile* tprof_recotruedx_truedx_yview = new TProfile("tprof_recotruedx_truedx_yview",
912  "Reco/True PCHit Path Length (Y View);True Path Length (cm);Reco/True Path Length",
913  tprof_bins,1,10);
914 
915  TProfile* tprof_dxres_truedx_xview = new TProfile("tprof_dxres_truedx_xview",
916  "Path Length Resolution (X View);True Path Length (cm);(Reco-True) / True Path Length",
917  tprof_bins,1,10);
918  TProfile* tprof_dxres_truedx_yview = new TProfile("tprof_dxres_truedx_yview",
919  "Path Length Resolution (Y View);True Path Length (cm);(Reco-True) / True Path Length",
920  tprof_bins,1,10);
921 
922  TProfile* tprof_dxres_w_xview = new TProfile("tprof_dxres_w_xview",
923  "Path Length Resolution (X View);W (cm);(Reco-True) / True Path Length",
924  tprof_bins,-wrange,wrange);
925  TProfile* tprof_dxres_w_yview = new TProfile("tprof_dxres_w_yview",
926  "Path Length Resolution (Y View);W (cm);(Reco-True) / True Path Length",
927  tprof_bins,-wrange,wrange);
928 
929 
930  TProfile* tprof_pecm_plane_xview = new TProfile("tprof_pecm_plane_xview",
931  "PE/cm vs Plane (X view);Plane;PE/cm",
932  PlanesPerView, 0, 2*PlanesPerView);
933  TProfile* tprof_pecm_plane_yview = new TProfile("tprof_pecm_plane_yview",
934  "PE/cm vs Plane (Y view);Plane;PE/cm",
935  PlanesPerView, 0, 2*PlanesPerView);
936 
937  TProfile* tprof_recodedx_dcosz_xview = new TProfile("tprof_recodedx_dcosz_xview",
938  "Reco dE/dx vs. Z Direction Cosine (X View);Z Direction Cosine;Reco dE/dx (MeV/cm)",
939  tprof_bins,-1,1);
940  TProfile* tprof_recodedx_dcosz_yview = new TProfile("tprof_recodedx_dcosz_yview",
941  "Reco dE/dx vs. Z Direction Cosine (Y View);Z Direction Cosine;Reco dE/dx (MeV/cm)",
942  tprof_bins,-1,1);
943  TProfile* tprof_truededx_dcosz_xview = new TProfile("tprof_truededx_dcosz_xview",
944  "True dE/dx vs. Z Direction Cosine (X View);Z Direction Cosine;True dE/dx (MeV/cm)",
945  tprof_bins,-1,1);
946  TProfile* tprof_truededx_dcosz_yview = new TProfile("tprof_truededx_dcosz_yview",
947  "True dE/dx vs. Z Direction Cosine (Y View);Z Direction Cosine;True dE/dx (MeV/cm)",
948  tprof_bins,-1,1);
949  TProfile* tprof_recodedx_dcosy_xview = new TProfile("tprof_recodedx_dcosy_xview",
950  "Reco dE/dx vs. Y Direction Cosine (X View);Y Direction Cosine;Reco dE/dx (MeV/cm)",
951  tprof_bins,-1,1);
952  TProfile* tprof_recodedx_dcosy_yview = new TProfile("tprof_recodedx_dcosy_yview",
953  "Reco dE/dx vs. Y Direction Cosine (Y View);Y Direction Cosine;Reco dE/dx (MeV/cm)",
954  tprof_bins,-1,1);
955  TProfile* tprof_truededx_dcosy_xview = new TProfile("tprof_truededx_dcosy_xview",
956  "True dE/dx vs. Y Direction Cosine (X View);Y Direction Cosine;True dE/dx (MeV/cm)",
957  tprof_bins,-1,1);
958  TProfile* tprof_truededx_dcosy_yview = new TProfile("tprof_truededx_dcosy_yview",
959  "True dE/dx vs. Y Direction Cosine (Y View);Y Direction Cosine;True dE/dx (MeV/cm)",
960  tprof_bins,-1,1);
961  TProfile* tprof_dx_dcosz_xview = new TProfile("tprof_dx_dcosz_xview",
962  "Reco Path Length vs. Z Direction Cos (X View);Z Direction Cosine;Reco Path Length (MeV/cm)",
963  tprof_bins,-1,1); // units wrong on X axis here, or title is wrong...
964  TProfile* tprof_dx_dcosz_yview = new TProfile("tprof_dx_dcosz_yview",
965  "Reco Path Length vs. Z Direction Cos (Y View);Z Direction Cosine;Reco Path Length (MeV/cm)",
966  tprof_bins,-1,1); // units wrong on X axis here, or title is wrong...
967  TProfile* tprof_dx_dcosy_xview = new TProfile("tprof_dx_dcosy_xview",
968  "Reco Path Length vs. Y Direction Cos (X View);Y Direction Cosine;Reco Path Length (MeV/cm)",
969  tprof_bins,-1,1); // units wrong on X axis here, or title is wrong...
970  TProfile* tprof_dx_dcosy_yview = new TProfile("tprof_dx_dcosy_yview",
971  "Reco Path Length vs. Y Direction Cos (Y View);Y Direction Cosine;Reco Path Length (MeV/cm)",
972  tprof_bins,-1,1); // units wrong on X axis here, or title is wrong...
973 
974 
975  /* 2D profiles
976 
977  tprof2_resp_cell-w_xview
978  tprof2_resp_w-cell_yview
979  tprof2_resp_cell-plane_*
980  tprof2_pecm_cell-plane_*
981 
982  */
983  // NOTE: left titles here alone!
984  TProfile2D *tprof2_resp_w_cell_xview = new TProfile2D("tprof2_resp_w-cell_xview",
985  "Mean Response (X view);X Cell;W (cm)",
986  nCells, 0, nCells,
987  tprof_bins, -wrange, wrange );
988  TProfile2D *tprof2_resp_cell_w_yview = new TProfile2D("tprof2_resp_cell-w_yview",
989  "Mean Response (Y yiew);W (cm);Y Cell",
990  tprof_bins, -wrange, wrange,
991  nCells, 0, nCells );
992 
993  TProfile2D *tprof2_resp_cell_plane_xview = new TProfile2D("tprof2_resp_cell-plane_xview",
994  "Detector Response PECorr/cm (X view);X Plane;Cell",
995  PlanesPerView, 0, 2*PlanesPerView,
996  nCells, 0, nCells);
997  TProfile2D *tprof2_resp_cell_plane_yview = new TProfile2D("tprof2_resp_cell-plane_yview",
998  "Detector Response PECorr/cm (Y view);Y Plane;Cell",
999  PlanesPerView, 0, 2*PlanesPerView,
1000  nCells, 0, nCells);
1001 
1002  TProfile2D *tprof2_pecm_cell_plane_xview = new TProfile2D("tprof2_pecm_cell-plane_xview",
1003  "Uncorrected Detector Response PE/cm (X view);X Plane;Cell",
1004  PlanesPerView, 0, 2*PlanesPerView,
1005  nCells, 0, nCells);
1006  TProfile2D *tprof2_pecm_cell_plane_yview = new TProfile2D("tprof2_pecm_cell-plane_yview",
1007  "Uncorrected Detector Response PE/cm (Y view);Y Plane;Cell",
1008  PlanesPerView, 0, 2*PlanesPerView,
1009  nCells, 0, nCells);
1010  TProfile2D *tprof2_pecm_cell_plane = new TProfile2D("tprof2_pecm_cell-plane",
1011  "Uncorrected Detector Response PE/cm;Plane;Cell",
1012  2*PlanesPerView, 0, 2*PlanesPerView,
1013  nCells, 0, nCells);
1014 
1015 
1016  // TH2D objects for uncalibrated cells
1017  TH2D *th2d_uncalib_cell_plane = new TH2D("th2d_uncalib_cell-plane",
1018  "Cells vs. Planes - Uncalibrated Cells;Plane;Cell",
1019  2*PlanesPerView, 0, 2*PlanesPerView,
1020  nCells, 0, nCells);
1021  TH2D *th2d_uncalib_cell_plane_xview = new TH2D("th2d_uncalib_cell-plane_xview",
1022  "Cells vs. Planes - Uncalibrated Cells (X view);Plane;Cell",
1023  2*PlanesPerView, 0, 2*PlanesPerView,
1024  nCells, 0, nCells);
1025  TH2D *th2d_uncalib_cell_plane_yview = new TH2D("th2d_uncalib_cell-plane_yview",
1026  "Cells vs. Planes - Uncalibrated Cells (Y view);Plane;Cell",
1027  2*PlanesPerView, 0, 2*PlanesPerView,
1028  nCells, 0, nCells);
1029 
1030 
1031  ////////////////////////////////////////////////////////////////////////
1032  // NOTE: New -- Now for the diblock versions... Do we want to do the range too?
1033  std::vector< std::vector< TProfile* > > tprof_resp_w_xview_fb_db; tprof_resp_w_xview_fb_db.resize(nFBbins+1);
1034  std::vector< std::vector< TProfile* > > tprof_resp_w_yview_fb_db; tprof_resp_w_yview_fb_db.resize(nFBbins+1);
1035  if( byDB ){
1036  for( unsigned int iFiber=0; iFiber < nFBbins+1; ++iFiber ){
1037  tprof_resp_w_xview_fb_db[iFiber].resize(nChan);
1038  tprof_resp_w_yview_fb_db[iFiber].resize(nChan);
1039  }
1040  }
1041  ////////////////////////////////////////////////////////////////////////
1042 
1043  /* fb-split
1044 
1045  N.B. -- all-FB plots kept in nFBbins index
1046 
1047  */
1048  // vs W
1049  std::vector< TProfile* > tprof_resp_w_xview_fb; tprof_resp_w_xview_fb.resize(nFBbins+1);
1050  std::vector< TProfile* > tprof_resp_w_yview_fb; tprof_resp_w_yview_fb.resize(nFBbins+1);
1051  std::vector< TProfile* > tprof_resp_w_xview_edgecells_fb; tprof_resp_w_xview_edgecells_fb.resize(nFBbins+1);
1052  std::vector< TProfile* > tprof_resp_w_yview_edgecells_fb; tprof_resp_w_yview_edgecells_fb.resize(nFBbins+1);
1053  std::vector< TH1D* > th1d_nhits_w_xview_fb; th1d_nhits_w_xview_fb.resize(nFBbins+1);
1054  std::vector< TH1D* > th1d_nhits_w_yview_fb; th1d_nhits_w_yview_fb.resize(nFBbins+1);
1055  std::vector< TH2D* > th2d_nhits_truededx_w_xview_fb; th2d_nhits_truededx_w_xview_fb.resize(nFBbins+1);
1056  std::vector< TH2D* > th2d_nhits_truededx_w_yview_fb; th2d_nhits_truededx_w_yview_fb.resize(nFBbins+1);
1057  std::vector< TH2D* > th2d_nhits_truededx_w_xview_nowindow_fb; th2d_nhits_truededx_w_xview_nowindow_fb.resize(nFBbins+1);
1058  std::vector< TH2D* > th2d_nhits_truededx_w_yview_nowindow_fb; th2d_nhits_truededx_w_yview_nowindow_fb.resize(nFBbins+1);
1059  std::vector< TProfile* > tprof_pe_w_xview; tprof_pe_w_xview.resize(nFBbins+1);
1060  std::vector< TProfile* > tprof_pe_w_yview; tprof_pe_w_yview.resize(nFBbins+1);
1061  std::vector< TProfile* > tprof_pecm_w_xview; tprof_pecm_w_xview.resize(nFBbins+1);
1062  std::vector< TProfile* > tprof_pecm_w_yview; tprof_pecm_w_yview.resize(nFBbins+1);
1063  std::vector< TProfile* > tprof_truededx_w_xview; tprof_truededx_w_xview.resize(nFBbins+1);
1064  std::vector< TProfile* > tprof_truededx_w_yview; tprof_truededx_w_yview.resize(nFBbins+1);
1065  std::vector< TProfile* > tprof_recodedx_w_xview; tprof_recodedx_w_xview.resize(nFBbins+1);
1066  std::vector< TProfile* > tprof_recodedx_w_yview; tprof_recodedx_w_yview.resize(nFBbins+1);
1067  std::vector< TProfile* > tprof_dx_w_xview; tprof_dx_w_xview.resize(nFBbins+1);
1068  std::vector< TProfile* > tprof_dx_w_yview; tprof_dx_w_yview.resize(nFBbins+1);
1069  std::vector< TProfile* > tprof_pecorrpe_w_xview_fb; tprof_pecorrpe_w_xview_fb.resize(nFBbins+1);
1070  std::vector< TProfile* > tprof_pecorrpe_w_yview_fb; tprof_pecorrpe_w_yview_fb.resize(nFBbins+1);
1071  std::vector< TProfile* > tprof_pecorrmev_w_xview_fb; tprof_pecorrmev_w_xview_fb.resize(nFBbins+1);
1072  std::vector< TProfile* > tprof_pecorrmev_w_yview_fb; tprof_pecorrmev_w_yview_fb.resize(nFBbins+1);
1073  // vs PE
1074  std::vector< TProfile* > tprof_resp_pe_xview_fb; tprof_resp_pe_xview_fb.resize(nFBbins+1);
1075  std::vector< TProfile* > tprof_resp_pe_yview_fb; tprof_resp_pe_yview_fb.resize(nFBbins+1);
1076  std::vector< TProfile* > tprof_pecorr_pe_xview_fb; tprof_pecorr_pe_xview_fb.resize(nFBbins+1);
1077  std::vector< TProfile* > tprof_pecorr_pe_yview_fb; tprof_pecorr_pe_yview_fb.resize(nFBbins+1);
1078  std::vector< TProfile* > tprof_pecorrmev_pe_xview_fb; tprof_pecorrmev_pe_xview_fb.resize(nFBbins+1);
1079  std::vector< TProfile* > tprof_pecorrmev_pe_yview_fb; tprof_pecorrmev_pe_yview_fb.resize(nFBbins+1);
1080  std::vector< TProfile* > tprof_pecorrpe_pe_xview_fb; tprof_pecorrpe_pe_xview_fb.resize(nFBbins+1);
1081  std::vector< TProfile* > tprof_pecorrpe_pe_yview_fb; tprof_pecorrpe_pe_yview_fb.resize(nFBbins+1);
1082  std::vector< TProfile* > tprof_truededx_pe_xview_fb; tprof_truededx_pe_xview_fb.resize(nFBbins+1);
1083  std::vector< TProfile* > tprof_truededx_pe_yview_fb; tprof_truededx_pe_yview_fb.resize(nFBbins+1);
1084  std::vector< TProfile* > tprof_recodedx_pe_xview_fb; tprof_recodedx_pe_xview_fb.resize(nFBbins+1);
1085  std::vector< TProfile* > tprof_recodedx_pe_yview_fb; tprof_recodedx_pe_yview_fb.resize(nFBbins+1);
1086  // vs PECorr
1087  std::vector< TProfile* > tprof_pecorrpe_pecorr_xview_fb; tprof_pecorrpe_pecorr_xview_fb.resize(nFBbins+1);
1088  std::vector< TProfile* > tprof_pecorrpe_pecorr_yview_fb; tprof_pecorrpe_pecorr_yview_fb.resize(nFBbins+1);
1089  std::vector< TProfile* > tprof_truededx_pecorr_xview_fb; tprof_truededx_pecorr_xview_fb.resize(nFBbins+1);
1090  std::vector< TProfile* > tprof_truededx_pecorr_yview_fb; tprof_truededx_pecorr_yview_fb.resize(nFBbins+1);
1091  std::vector< TProfile* > tprof_recodedx_pecorr_xview_fb; tprof_recodedx_pecorr_xview_fb.resize(nFBbins+1);
1092  std::vector< TProfile* > tprof_recodedx_pecorr_yview_fb; tprof_recodedx_pecorr_yview_fb.resize(nFBbins+1);
1093  // vs true dE
1094  std::vector< TProfile* > tprof_pecorrmev_dE_xview_fb; tprof_pecorrmev_dE_xview_fb.resize(nFBbins+1);
1095  std::vector< TProfile* > tprof_pecorrmev_dE_yview_fb; tprof_pecorrmev_dE_yview_fb.resize(nFBbins+1);
1096  // vs x (dist from end of track)
1097  std::vector< TProfile* > tprof_truededx_x_xview; tprof_truededx_x_xview.resize(nFBbins+1);
1098  std::vector< TProfile* > tprof_truededx_x_yview; tprof_truededx_x_yview.resize(nFBbins+1);
1099  std::vector< TProfile* > tprof_recodedx_x_xview; tprof_recodedx_x_xview.resize(nFBbins+1);
1100  std::vector< TProfile* > tprof_recodedx_x_yview; tprof_recodedx_x_yview.resize(nFBbins+1);
1101  // vs pathlength
1102  std::vector< TProfile* > tprof_truededx_dx_xview; tprof_truededx_dx_xview.resize(nFBbins+1);
1103  std::vector< TProfile* > tprof_truededx_dx_yview; tprof_truededx_dx_yview.resize(nFBbins+1);
1104  std::vector< TProfile* > tprof_recodedx_dx_xview; tprof_recodedx_dx_xview.resize(nFBbins+1);
1105  std::vector< TProfile* > tprof_recodedx_dx_yview; tprof_recodedx_dx_yview.resize(nFBbins+1);
1106  // vs cells per plane
1107  std::vector< TProfile* > tprof_truededx_cellpplane_xview; tprof_truededx_cellpplane_xview.resize(nFBbins+1);
1108  std::vector< TProfile* > tprof_truededx_cellpplane_yview; tprof_truededx_cellpplane_yview.resize(nFBbins+1);
1109  std::vector< TProfile* > tprof_recodedx_cellpplane_xview; tprof_recodedx_cellpplane_xview.resize(nFBbins+1);
1110  std::vector< TProfile* > tprof_recodedx_cellpplane_yview; tprof_recodedx_cellpplane_yview.resize(nFBbins+1);
1111  // vs FB
1112  std::vector< TH1D* > th1d_nhits_resp_xview_fb; th1d_nhits_resp_xview_fb.resize(nFBbins+1);
1113  std::vector< TH1D* > th1d_nhits_resp_yview_fb; th1d_nhits_resp_yview_fb.resize(nFBbins+1);
1114  std::vector< TH1D* > th1d_nhits_truededx_xview_fb; th1d_nhits_truededx_xview_fb.resize(nFBbins+1);
1115  std::vector< TH1D* > th1d_nhits_truededx_yview_fb; th1d_nhits_truededx_yview_fb.resize(nFBbins+1);
1116  std::vector< TH1D* > th1d_nhits_truedetruedx_xview_fb; th1d_nhits_truedetruedx_xview_fb.resize(nFBbins+1);
1117  std::vector< TH1D* > th1d_nhits_truedetruedx_yview_fb; th1d_nhits_truedetruedx_yview_fb.resize(nFBbins+1);
1118  std::vector< TH1D* > th1d_nhits_recodedx_xview_fb; th1d_nhits_recodedx_xview_fb.resize(nFBbins+1);
1119  std::vector< TH1D* > th1d_nhits_recodedx_yview_fb; th1d_nhits_recodedx_yview_fb.resize(nFBbins+1);
1120  std::vector< TH1D* > th1d_nhits_pe_xview_fb; th1d_nhits_pe_xview_fb.resize(nFBbins+1);
1121  std::vector< TH1D* > th1d_nhits_pe_yview_fb; th1d_nhits_pe_yview_fb.resize(nFBbins+1);
1122  std::vector< TH1D* > th1d_nhits_pecm_xview_fb; th1d_nhits_pecm_xview_fb.resize(nFBbins+1);
1123  std::vector< TH1D* > th1d_nhits_pecm_yview_fb; th1d_nhits_pecm_yview_fb.resize(nFBbins+1);
1124  std::vector< TH2D* > th2d_nhits_pecm_w_xview_fb; th2d_nhits_pecm_w_xview_fb.resize(nFBbins+1);
1125  std::vector< TH2D* > th2d_nhits_pecm_w_yview_fb; th2d_nhits_pecm_w_yview_fb.resize(nFBbins+1);
1126  std::string name, fbtag, titleC, titleX, titleY;
1127  //
1128  // plotname convention: "plottype_fillweight_axis(-axis)_specifiers"
1129  //
1130  // Specifiers will usually be view and fb bin, so make name string:
1131  // "plottype_fillweight_axis(-axis)"
1132  // and handle view/fb with Form()
1133  //
1134  // UPDATED CONVENTION!!!
1135  // 1. In the plot names, label the specifier for FB as follows:
1136  // a. fb < nFBbins -> "-fbsplit{fb%i}"
1137  // b. else "" (leave this alone so that --I THINK-- Tyler's calc_E_scale macro will still pick up histograms)
1138  // 2. In the titles, use a Form() as well, with titleC, titleX, and titleY as strings
1139  // a. convention: "Y axis vs. X axis (X/Y view)"
1140  //
1141  for(unsigned int i=0; i<nFBbins+1; i++){
1142  if( i < nFBbins ) fbtag = Form("-fbsplit{fb%i}",i);
1143  else fbtag="";
1144 
1145  name = "tprof_resp_w_";
1146  titleC = "Mean Calibrated Response vs. W";
1147  titleX = "W (cm)";
1148  titleY = "Response (PECorr/cm)";
1149  tprof_resp_w_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1150  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1151  tprof_bins,-wrange,wrange);
1152  tprof_resp_w_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1153  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1154  tprof_bins,-wrange,wrange);
1155 
1156  ////////////////////////////////////////////////////////////////////////////////////////
1157  // NOTE: New -- Now for the diblock versions... Do we want range split too?
1158  if( byDB ){
1159  for( unsigned int thisDB=0; thisDB<nChan; ++thisDB ){
1160  std::string dbtag;
1161  std::string othertag = Form("db%u_",thisDB+1);
1162  if( i < nFBbins )
1163  dbtag = Form("-fbsplit{fb%i}",i);
1164  else
1165  dbtag = "";
1166  tprof_resp_w_xview_fb_db[i][thisDB] = new TProfile(Form("%s%sxview%s",othertag.c_str(),name.c_str(),dbtag.c_str()),
1167  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1168  tprof_bins,-wrange,wrange);
1169  tprof_resp_w_yview_fb_db[i][thisDB] = new TProfile(Form("%s%syview%s",othertag.c_str(),name.c_str(),dbtag.c_str()),
1170  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1171  tprof_bins,-wrange,wrange);
1172  }
1173  }
1174  ////////////////////////////////////////////////////////////////////////////////////////
1175 
1176 
1177  name = "tprof_resp_w_edgecells-";
1178  titleC = "Mean Calibrated Response in Edge Cells vs. W";
1179  titleX = "W (cm)";
1180  titleY = "Response (PECorr/cm)";
1181  tprof_resp_w_xview_edgecells_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1182  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1183  tprof_bins,-wrange,wrange);
1184  tprof_resp_w_yview_edgecells_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1185  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1186  tprof_bins,-wrange,wrange);
1187 
1188  name = "th1d_nhits_w_";
1189  titleC = "Number of Hits vs. W";
1190  titleX = "W (cm)";
1191  titleY = "# Hits";
1192  th1d_nhits_w_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1193  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1194  tprof_bins,-wrange,wrange);
1195  th1d_nhits_w_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1196  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1197  tprof_bins,-wrange,wrange);
1198 
1199  name = "th2d_nhits_truededx-w_";
1200  titleC = "True Response vs. W (nhits)"; // is this the best way to address TH2Ds in the convention stated? Going to follow this.
1201  titleX = "W (cm)";
1202  titleY = "True Response (MeV/cm)";
1203  th2d_nhits_truededx_w_xview_fb[i] = new TH2D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1204  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1205  tprof_bins,-wrange,wrange,
1206  tprof_bins,0,10);
1207  th2d_nhits_truededx_w_yview_fb[i] = new TH2D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1208  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1209  tprof_bins,-wrange,wrange,
1210  tprof_bins,0,10);
1211 
1212  name = "th2d_nhits_truededx-w_nowindow-";
1213  titleC = "True Response vs. W (hits not in track window)";
1214  titleX = "W (cm)";
1215  titleY = "True Response (MeV/cm)";
1216  th2d_nhits_truededx_w_xview_nowindow_fb[i] = new TH2D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1217  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1218  tprof_bins,-wrange,wrange,
1219  tprof_bins,0,10);
1220  th2d_nhits_truededx_w_yview_nowindow_fb[i] = new TH2D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1221  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1222  tprof_bins,-wrange,wrange,
1223  tprof_bins,0,10);
1224 
1225  name = "tprof_pe_w_";
1226  titleC = "PE vs. W (hits not in track window)";
1227  titleX = "W (cm)";
1228  titleY = "PE";
1229  tprof_pe_w_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1230  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1231  tprof_bins,-wrange,wrange);
1232  tprof_pe_w_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1233  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1234  tprof_bins,-wrange,wrange);
1235 
1236 
1237  name = "tprof_pecm_w_";
1238  titleC = "Uncorrected Response vs. W";
1239  titleX = "W (cm)";
1240  titleY = "Uncorrected Response (PE/cm)";
1241  tprof_pecm_w_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1242  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1243  tprof_bins,-wrange,wrange);
1244  tprof_pecm_w_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1245  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1246  tprof_bins,-wrange,wrange);
1247 
1248  name = "tprof_truededx_w_";
1249  titleC = "True dE/dx vs. W";
1250  titleX = "W (cm)";
1251  titleY = "True dE/dx (MeV/cm)";
1252  tprof_truededx_w_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1253  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1254  tprof_bins,-wrange,wrange);
1255  tprof_truededx_w_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1256  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1257  tprof_bins,-wrange,wrange);
1258 
1259  name = "tprof_recodedx_w_";
1260  titleC = "Reco dE/dx vs. W";
1261  titleX = "W (cm)";
1262  titleY = "Reco dE/dx (MeV/cm)";
1263  tprof_recodedx_w_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1264  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1265  tprof_bins,-wrange,wrange);
1266  tprof_recodedx_w_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1267  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1268  tprof_bins,-wrange,wrange);
1269 
1270  name = "tprof_dx_w_";
1271  titleC = "Path Length vs. W";
1272  titleX = "W (cm)";
1273  titleY = "Path Length (cm)";
1274  tprof_dx_w_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1275  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1276  tprof_bins,-wrange,wrange);
1277  tprof_dx_w_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1278  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1279  tprof_bins,-wrange,wrange);
1280 
1281  // Check this is indeed meant to be vs. W
1282  name = "tprof_pecorrpe_w_";
1283  titleC = "PECorr/PE vs. W";
1284  titleX = "W (cm)";
1285  titleY = "PECorr/PE";
1286  tprof_pecorrpe_w_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1287  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1288  tprof_bins,-wrange,wrange);
1289  tprof_pecorrpe_w_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1290  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1291  tprof_bins,-wrange,wrange);
1292 
1293  name = "tprof_pecorrmev_w_";
1294  titleC = "PECorr/MeV vs. W";
1295  titleX = "W (cm)";
1296  titleY = "PECorr/MeV";
1297  tprof_pecorrmev_w_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1298  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1299  tprof_bins,-wrange,wrange);
1300  tprof_pecorrmev_w_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1301  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1302  tprof_bins,-wrange,wrange);
1303 
1304  name = "tprof_resp_pe_";
1305  titleC = "Corrected Response (PECorr/cm) vs. PE";
1306  titleX = "PE";
1307  titleY = "Corrected Response (PECorr/cm)";
1308  tprof_resp_pe_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1309  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1310  /*"Response vs PE (X view);PE;PECorr/cm",*/
1311  tprof_bins,pe_x_min,pe_x_max);
1312  tprof_resp_pe_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1313  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1314  /*"Response vs PE (Y view);PE;PECorr/cm",*/
1315  tprof_bins,pe_y_min,pe_y_max);
1316 
1317  name = "tprof_pecorr_pe_";
1318  titleC = "PECorr vs. PE";
1319  titleX = "PE";
1320  titleY = "PECorr";
1321  tprof_pecorr_pe_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1322  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1323  tprof_bins,pe_x_min,pe_x_max);
1324  tprof_pecorr_pe_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1325  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1326  tprof_bins,pe_y_min,pe_y_max);
1327 
1328  name = "tprof_pecorrmev_pe_";
1329  titleC = "PECorr/MeV vs. PE";
1330  titleX = "PE";
1331  titleY = "PECorr/MeV";
1332  tprof_pecorrmev_pe_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1333  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1334  tprof_bins,pe_x_min,pe_x_max);
1335  tprof_pecorrmev_pe_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1336  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1337  tprof_bins,pe_y_min,pe_y_max);
1338 
1339  name = "tprof_pecorrpe_pe_";
1340  titleC = "PECorr/PE vs. PE";
1341  titleX = "PE";
1342  titleY = "PECorr/PE";
1343  tprof_pecorrpe_pe_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1344  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1345  tprof_bins,pe_x_min,pe_x_max);
1346  tprof_pecorrpe_pe_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1347  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1348  tprof_bins,pe_y_min,pe_y_max);
1349 
1350  name = "tprof_truededx_pe_";
1351  titleC = "True dE/dx vs. PE";
1352  titleX = "PE";
1353  titleY = "True dE/dx (MeV/cm)";
1354  tprof_truededx_pe_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1355  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1356  tprof_bins,pe_y_min,pe_y_max);
1357  tprof_truededx_pe_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1358  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1359  tprof_bins,pe_y_min,pe_y_max);
1360 
1361  name = "tprof_recodedx_pe_";
1362  titleC = "Reco dE/dx vs. PE";
1363  titleX = "PE";
1364  titleY = "Reco dE/dx (MeV/cm)";
1365  tprof_recodedx_pe_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1366  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1367  tprof_bins,pe_y_min,pe_y_max);
1368  tprof_recodedx_pe_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1369  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1370  tprof_bins,pe_y_min,pe_y_max);
1371 
1372  name = "tprof_truededx_pecorr_";
1373  titleC = "True dE/dx vs. PECorr";
1374  titleX = "PECorr";
1375  titleY = "True dE/dx (MeV/cm)";
1376  tprof_truededx_pecorr_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1377  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1378  tprof_bins,pe_y_min,pe_y_max);
1379  tprof_truededx_pecorr_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1380  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1381  tprof_bins,pe_y_min,pe_y_max);
1382 
1383  name = "tprof_recodedx_pecorr_";
1384  titleC = "Reco dE/dx vs. PECorr";
1385  titleX = "PECorr";
1386  titleY = "Reco dE/dx (MeV/cm)";
1387  tprof_recodedx_pecorr_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1388  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1389  tprof_bins,pe_y_min,pe_y_max);
1390  tprof_recodedx_pecorr_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1391  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1392  tprof_bins,pe_y_min,pe_y_max);
1393 
1394  name = "tprof_pecorrpe_pecorr_";
1395  titleC = "PECorr/PE vs. PECorr";
1396  titleX = "PECorr";
1397  titleY = "PECorr/PE";
1398  tprof_pecorrpe_pecorr_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1399  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1400  tprof_bins,pe_y_min,pe_y_max);
1401  tprof_pecorrpe_pecorr_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1402  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1403  tprof_bins,pe_y_min,pe_y_max);
1404 
1405 
1406  name = "tprof_truededx_x_";
1407  titleC = "True dE/dx vs. distance to end of track";
1408  titleX = "Distance to End of Track";
1409  titleY = "True dE/dx (MeV/cm)";
1410  tprof_truededx_x_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1411  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1412  tprof_bins,0,track_end_max);
1413  tprof_truededx_x_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1414  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1415  tprof_bins,0,track_end_max);
1416 
1417  name = "tprof_recodedx_x_";
1418  titleC = "Reco dE/dx vs. distance to end of track";
1419  titleX = "Distance to End of Track";
1420  titleY = "Reco dE/dx (MeV/cm)";
1421  tprof_recodedx_x_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1422  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1423  tprof_bins,0,track_end_max);
1424  tprof_recodedx_x_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1425  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1426  tprof_bins,0,track_end_max);
1427 
1428  name = "tprof_truededx_dx_";
1429  titleC = "True dE/dx vs. Reco Path Length";
1430  titleX = "Reco Path Length (cm)";
1431  titleY = "True dE/dx (MeV/cm)";
1432  tprof_truededx_dx_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1433  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1434  dxbins,dxmin,dxmax);
1435  tprof_truededx_dx_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1436  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1437  dxbins,dxmin,dxmax);
1438 
1439  name = "tprof_recodedx_dx_";
1440  titleC = "Reco dE/dx vs. Reco Path Length";
1441  titleX = "Reco Path Length (cm)";
1442  titleY = "Reco dE/dx (MeV/cm)";
1443  tprof_recodedx_dx_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1444  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1445  dxbins,dxmin,dxmax);
1446  tprof_recodedx_dx_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1447  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1448  dxbins,dxmin,dxmax);
1449 
1450  name = "tprof_truededx_cellpplane";
1451  titleC = "True dE/dx vs. Number Cells Per Plane";
1452  titleX = "# of Cells Per Plane";
1453  titleY = "True dE/dx (MeV/cm)";
1454  tprof_truededx_cellpplane_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1455  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1456  10,0,10);
1457  tprof_truededx_cellpplane_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1458  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1459  10,0,10);
1460 
1461  name = "tprof_recodedx_cellpplane";
1462  titleC = "Reco dE/dx vs. Number Cells Per Plane";
1463  titleX = "# of Cells Per Plane";
1464  titleY = "Reco dE/dx (MeV/cm)";
1465  tprof_recodedx_cellpplane_xview[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1466  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1467  10,0,10);
1468  tprof_recodedx_cellpplane_yview[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1469  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1470  10,0,10);
1471 
1472  name = "tprof_pecorrmev_dE_";
1473  titleC = "PECorr/MeV vs. True Energy Deposit";
1474  titleX = "True Energy Deposit (MeV)";
1475  titleY = "PECorr/MeV";
1476  tprof_pecorrmev_dE_xview_fb[i] = new TProfile(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1477  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1478  tprof_bins,3,25);
1479  tprof_pecorrmev_dE_yview_fb[i] = new TProfile(Form("%syview%s",name.c_str(),fbtag.c_str()),
1480  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1481  tprof_bins,3,25);
1482 
1483 
1484  name = "th1d_nhits_resp_";
1485  titleC = "N Hits vs. Reconstructed Response"; // Proper title??
1486  titleX = "Reconstructed Response (PECorr/cm)";
1487  titleY = "N hits"; // Is this the proper y axis??
1488  th1d_nhits_resp_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1489  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1490  /*"Reconstructed Response (X View);PECorr/cm;",*/
1491  tprof_bins,0,100);
1492  th1d_nhits_resp_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1493  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1494  /*"Reconstructed Response (Y View);PECorr/cm;",*/
1495  tprof_bins,0,100);
1496 
1497  name = "th1d_nhits_truededx_";
1498  titleC = "N Hits vs. True Response"; // Proper title and axis??
1499  titleX = "True Response (MeV/cm)";
1500  titleY = "N hits";
1501  th1d_nhits_truededx_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1502  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1503  /*"True Response (X View);MeV/cm;",*/
1504  tprof_bins,0,10);
1505  th1d_nhits_truededx_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1506  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1507  /*"True Response (Y View);MeV/cm;",*/
1508  tprof_bins,0,10);
1509 
1510  name = "th1d_nhits_truedetruedx_";
1511  titleC = "N Hits vs. True Response"; // proper title?
1512  titleX = "True Response (true MeV / true cm)";
1513  titleY = "N hits";
1514  th1d_nhits_truedetruedx_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1515  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1516  /*"True Response (X View);true MeV / true cm;",*/
1517  tprof_bins,0,4);
1518  th1d_nhits_truedetruedx_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1519  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1520  /*"True Response (Y View);true MeV / true cm;",*/
1521  tprof_bins,0,4);
1522 
1523  name = "th1d_nhits_recodedx_";
1524  titleC = "N Hits vs. Reco dE/dx"; // NOTE: response and dE/dx interchanged fairly liberally in these titles. Should we make uniform?
1525  titleX = "Reco dE/dx (MeV/cm)";
1526  titleY = "N hits";
1527  th1d_nhits_recodedx_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1528  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1529  /*"X View;Reco dE/dx (MeV/cm);",*/
1530  tprof_bins,0,10);
1531  th1d_nhits_recodedx_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1532  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1533  /*"Y View;Reco dE/dx (MeV/cm);",*/
1534  tprof_bins,0,10);
1535 
1536  name = "th1d_nhits_pe_";
1537  titleC = "N Hits vs. PE"; // again, proper title?
1538  titleX = "PE";
1539  titleY = "N hits";
1540  th1d_nhits_pe_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1541  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1542  /*"X View;PE;",*/
1543  tprof_bins,0,450);
1544  th1d_nhits_pe_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1545  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1546  /*"Y View;PE;",*/
1547  tprof_bins,0,450);
1548 
1549  name = "th1d_nhits_pecm_";
1550  titleC = "N Hits vs. Uncorrected Response"; // same with response and PE/cm - do we want to be uniform here?
1551  titleX = "Uncorrected Response (PE/cm)";
1552  titleY = "N hits";
1553  th1d_nhits_pecm_xview_fb[i] = new TH1D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1554  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1555  /*"X View;PE/cm;",*/
1556  tprof_bins,0,100);
1557  th1d_nhits_pecm_yview_fb[i] = new TH1D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1558  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1559  /*"Y View;PE/cm;",*/
1560  tprof_bins,0,100);
1561 
1562  name = "th2d_nhits_pecm-w_";
1563  titleC = "W vs. Uncorrected Response (N hits)";
1564  titleX = "Uncorrected Response (PE/cm)";
1565  titleY = "W (cm)";
1566  th2d_nhits_pecm_w_xview_fb[i] = new TH2D(Form("%sxview%s",name.c_str(),fbtag.c_str()),
1567  Form("%s (X view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1568  tprof_bins, -wrange, wrange,
1569  tprof_bins, pe_x_min, pe_x_max );
1570  th2d_nhits_pecm_w_yview_fb[i] = new TH2D(Form("%syview%s",name.c_str(),fbtag.c_str()),
1571  Form("%s (Y view);%s;%s",titleC.c_str(),titleX.c_str(),titleY.c_str()),
1572  tprof_bins, -wrange, wrange,
1573  tprof_bins, pe_y_min, pe_y_max );
1574  }
1575 
1576  std::cout << "...plots declared" << std::endl;
1577 
1578  const unsigned int nWslices = 8;
1579  std::vector<float> minW; minW.resize(nWslices);
1580  std::vector<float> maxW; maxW.resize(nWslices);
1581  if(det=="nd"){
1582  minW[0] = -201;
1583  minW[1] = -150;
1584  minW[2] = -100;
1585  minW[3] = -50;
1586  minW[4] = 0;
1587  minW[5] = 50;
1588  minW[6] = 100;
1589  minW[7] = 150;
1590  maxW[0] = -150;
1591  maxW[1] = -100;
1592  maxW[2] = -50;
1593  maxW[3] = 0;
1594  maxW[4] = 50;
1595  maxW[5] = 100;
1596  maxW[6] = 150;
1597  maxW[7] = 200;
1598  }
1599  if(det=="fd"){
1600  minW[0] = -800;
1601  minW[1] = -600;
1602  minW[2] = -400;
1603  minW[3] = -200;
1604  minW[4] = 0;
1605  minW[5] = 200;
1606  minW[6] = 400;
1607  minW[7] = 600;
1608  maxW[0] = -600;
1609  maxW[1] = -400;
1610  maxW[2] = -200;
1611  maxW[3] = 0;
1612  maxW[4] = 200;
1613  maxW[5] = 400;
1614  maxW[6] = 600;
1615  maxW[7] = 800;
1616  }
1617 
1618  std::vector< TProfile* > tprof_pe_truedE_xview_w; tprof_pe_truedE_xview_w.resize(nWslices);
1619  std::vector< TProfile* > tprof_pe_truedE_yview_w; tprof_pe_truedE_yview_w.resize(nWslices);
1620  for(unsigned int i=0; i<nWslices; i++){
1621 
1622  tprof_pe_truedE_xview_w[i] = new TProfile(Form("tprof_pe_truedE_xview_w.%i.%i",(int)(minW[i]),(int)(maxW[i])),
1623  Form("Threshold Effect, W = %i to %i;PE;True Deposition (MeV)",
1624  (int)(minW[i]),(int)(maxW[i])),
1625  tprof_bins,0,40);
1626  tprof_pe_truedE_yview_w[i] = new TProfile(Form("tprof_pe_truedE_yview_w.%i.%i",(int)(minW[i]),(int)(maxW[i])),
1627  Form("Threshold Effect, W = %i to %i;PE;True Deposition (MeV)",
1628  (int)(minW[i]),(int)(maxW[i])),
1629  tprof_bins,0,40);
1630  }
1631 
1632  ///////////////////////////////////////////////////////////////////
1633  ///////////////////////////////////////////////////////////////////
1634  ////////////////// //////////////////
1635  ////////////////// Loop through TTree entries //////////////////
1636  ////////////////// //////////////////
1637  ///////////////////////////////////////////////////////////////////
1638  ///////////////////////////////////////////////////////////////////
1639  unsigned int n = t->GetEntries();
1640  unsigned int n_zero_true_path=0;
1641  unsigned int n_zero_path=0;
1642  unsigned int n_zero_PEc=0;
1643  unsigned int n_neg_PEc=0;
1644  unsigned int n_zero_trueE=0;
1645  unsigned int n_track_window_hits=0;
1646  unsigned int n_nonmuon_hits=0;
1647  //if (n > 10e6) n = 10e6;
1648  if (maxEntries > 0 && n > maxEntries) n = maxEntries;
1649  std::cout << "~~~~ looping through "
1650  << n << " TTree entries ~~~~"
1651  << std::endl;
1652  for( unsigned int e=0; e<n; e++ ){
1653  t->GetEntry(e);
1654 
1655  if(tree_i==1||tree_i==3){
1656  if(view==0) w=ypos;
1657  if(view==1) w=xpos;
1658  dcosz = zDir / std::sqrt(xDir*xDir + yDir*yDir + zDir*zDir);
1659  dcosy = yDir / std::sqrt(xDir*xDir + yDir*yDir + zDir*zDir);
1660  }
1661 
1662  if(tree_i==0){
1663  cell = (int)(cell_fl);
1664  plane = (int)(plane_fl);
1665  fbBin = (int)(fbBin_fl);
1666  }
1667 
1668  // gev-->mev
1669  if(tree_i!=0) trueE = 1000*trueE;
1670  float recoMeV = 1000*recoGeV;
1671  // this sample is negative??
1672  //if (sample.find("traj")!=std::string::npos) recoMeV = -1*recoMeV;
1673 
1674  // FD MC real conditions is different??
1675  //recoMeV = -1*recoGeV;
1676 
1677  if(view==0){
1678  th2d_nhits_dE_pe_xview->Fill(PE,trueE);
1679  }else{
1680  th2d_nhits_dE_pe_yview->Fill(PE,trueE);
1681  }
1682 
1683  if(e%1000000==0) std::cout << e << " / " << n << std::endl;
1684 
1685  // be careful -- inside or outside of track window cut
1686  //if( path <= 0 ) std::cout << "path = " << path << std::endl;
1687  //if( PE <= 0 ) std::cout << "PE = " << PE << std::endl;
1688  //if( true_path <= 0 && isMC ) n_zero_true_path++;
1689 
1690  int bin=fbBin;
1691  int all=nFBbins;
1692  if(bin<0||nFBbins-1<(unsigned int)bin) bin=0;
1693 
1694 
1695  // First fill things regardless of track window
1696  // -For very few things. better to just have same cut on all plots for one run
1697  // so it doesn't have to be constantly referenced which has what cut
1698  // -Just fill "nowindow" plots and track window plot
1699  if(view==0){
1700  th2d_nhits_truededx_w_xview_nowindow_fb [bin]->Fill( w, trueE/path );
1701  th2d_nhits_truededx_w_xview_nowindow_fb [all]->Fill( w, trueE/path );
1702  tprof_truededx_x_xview [bin]->Fill( x, trueE/path );
1703  tprof_truededx_x_xview [all]->Fill( x, trueE/path );
1704  tprof_recodedx_x_xview [bin]->Fill( x, recoMeV/path );
1705  tprof_recodedx_x_xview [all]->Fill( x, recoMeV/path );
1706  } else {
1707  th2d_nhits_truededx_w_yview_nowindow_fb [bin]->Fill( w, trueE/path );
1708  th2d_nhits_truededx_w_yview_nowindow_fb [all]->Fill( w, trueE/path );
1709  tprof_truededx_x_yview [bin]->Fill( x, trueE/path );
1710  tprof_truededx_x_yview [all]->Fill( x, trueE/path );
1711  tprof_recodedx_x_yview [bin]->Fill( x, recoMeV/path );
1712  tprof_recodedx_x_yview [all]->Fill( x, recoMeV/path );
1713  }
1714 
1715  bool inTrackWindow = cut_Track_Window(x);
1716  if( PECorr < 0. && th2d_uncalib_cell_plane->GetBinContent(plane+1,cell+1)==0. ){
1717  // Fill this cell to mark it as uncalibrated for validation purposes
1718  th2d_uncalib_cell_plane->Fill(plane,cell);
1719  if( view==0 ) th2d_uncalib_cell_plane_xview->Fill(plane,cell);
1720  if( view==1 ) th2d_uncalib_cell_plane_yview->Fill(plane,cell);
1721  }
1722  // Check between 100 and 200cm from end of track and Calibrator returns -1 for hits in uncalibrated cells and check PE>0
1723  // NOTE: THIS IS WHERE THE PLOTS USED NOW FOR CALIBRATION (AS OF ANA2017 and ANA2018, FLAT W HISTOGRAMS) ARE FILLED!!!!
1724  if( inTrackWindow && PECorr > 0. && cut_PE(PE) ){
1725  n_track_window_hits++;
1726 
1727  if(std::abs(pdg)!=13) n_nonmuon_hits++;
1728 
1729  if(PECorr==0.){
1730  n_zero_PEc++;
1731  }
1732  if(PECorr<=0.){
1733  n_neg_PEc++;
1734  }
1735  if(trueE<=0.){
1736  n_zero_trueE++;
1737  }
1738  if(path<=0.){
1739  n_zero_path++;
1740  }
1741 
1742  // flat W cuts for Absolute Calibration
1743  if( (view==0 && flatW_x_min < w && w < flatW_x_max) ||
1744  (view==1 && flatW_y_min < w && w < flatW_y_max) ){
1745  th1d_nhits_resp_flatW ->Fill(PECorr/path);
1746  th1d_nhits_truededx_flatW->Fill(trueE/path);
1747  if( nDivisions>0 ){
1748  th1d_nhits_resp_flatW_range[int(float(run-loRun)/float(runsPerDiv))] ->Fill(PECorr/path);
1749  }
1750  if( byDB ){
1751  th1d_nhits_resp_flatW_diblock[diblock-1] ->Fill(PECorr/path);
1752  th1d_nhits_truededx_flatW_diblock[diblock-1]->Fill(trueE/path);
1753  if( nDivisions>0 ){
1754  th1d_nhits_resp_flatW_range_diblock[int(float(run-loRun)/float(runsPerDiv))][diblock-1] ->Fill(PECorr/path);
1755  }
1756  }
1757  }
1758  if(view==0 && flatW_x_min < w && w < flatW_x_max){
1759  th1d_nhits_resp_flatW_xview ->Fill(PECorr/path);
1760  th1d_nhits_truededx_flatW_xview->Fill(trueE/path);
1761  th1d_nhits_recodedx_flatW_xview->Fill(recoMeV/path);
1762  if( tree_i==0 ) th1d_nxyhits_flatW_xview->Fill(nXYHits,1./double(nXYHits));
1763  if( nDivisions>0 ){
1764  th1d_nhits_resp_flatW_xview_range[int(float(run-loRun)/float(runsPerDiv))] ->Fill(PECorr/path);
1765  }
1766  if( byDB ){
1767  th1d_nhits_resp_flatW_xview_diblock[diblock-1] ->Fill(PECorr/path);
1768  th1d_nhits_truededx_flatW_xview_diblock[diblock-1]->Fill(trueE/path);
1769  if( nDivisions>0 ){
1770  th1d_nhits_resp_flatW_xview_range_diblock[int(float(run-loRun)/float(runsPerDiv))][diblock-1] ->Fill(PECorr/path);
1771  }
1772  }
1773  }
1774  if(view==1 && flatW_y_min < w && w < flatW_y_max){
1775  th1d_nhits_resp_flatW_yview ->Fill(PECorr/path);
1776  th1d_nhits_truededx_flatW_yview->Fill(trueE/path);
1777  th1d_nhits_recodedx_flatW_yview->Fill(recoMeV/path);
1778  if( tree_i==0 ) th1d_nxyhits_flatW_yview->Fill(nXYHits,1./double(nXYHits));
1779  if( nDivisions>0 ){
1780  th1d_nhits_resp_flatW_yview_range[int(float(run-loRun)/float(runsPerDiv))] ->Fill(PECorr/path);
1781  }
1782  if( byDB ){
1783  th1d_nhits_resp_flatW_yview_diblock[diblock-1] ->Fill(PECorr/path);
1784  th1d_nhits_truededx_flatW_yview_diblock[diblock-1]->Fill(trueE/path);
1785  if( nDivisions>0 ){
1786  th1d_nhits_resp_flatW_yview_range_diblock[int(float(run-loRun)/float(runsPerDiv))][diblock-1] ->Fill(PECorr/path);
1787  }
1788  }
1789  }
1790 
1791  // Fill the by-time plots here for reference, if tree_i==0
1792  // NEW!!!!
1793  if( tree_i==0 ){
1794  if(view==0 && flatW_x_min < w && w < flatW_x_max){
1795  tprof_resp_time_flatW_xview->Fill(EventTimeHigh,PECorr/path);
1796  tprof_pecm_time_flatW_xview->Fill(EventTimeHigh,PE/path);
1797  tprof_pe_time_flatW_xview->Fill(EventTimeHigh,PE);
1798  if( byDB ){
1799  tprof_pecm_time_flatW_xview_diblock[diblock-1]->Fill(EventTimeHigh,PE/path);
1800  }
1801  // W-breakdowns
1802  if( w>=flatW_x_min && w<flatW_x_min+flatW_x_selRgn ){
1803  //std::cout << "LO W = " << w << std::endl;
1804  tprof_resp_time_flatW_xview_loW->Fill(EventTimeHigh,PECorr/path);
1805  tprof_pecm_time_flatW_xview_loW->Fill(EventTimeHigh,PE/path);
1806  }
1807  else if( w>=flatW_x_min+flatW_x_selRgn && w<flatW_x_min+(2.*flatW_x_selRgn) ){
1808  //std::cout << "MD W = " << w << std::endl;
1809  tprof_resp_time_flatW_xview_mdW->Fill(EventTimeHigh,PECorr/path);
1810  tprof_pecm_time_flatW_xview_mdW->Fill(EventTimeHigh,PE/path);
1811  }
1812  else if( w>=flatW_x_min+(2.*flatW_x_selRgn) && w<flatW_x_min+(3.*flatW_x_selRgn) ){
1813  //std::cout << "HI W = " << w << std::endl;
1814  tprof_resp_time_flatW_xview_hiW->Fill(EventTimeHigh,PECorr/path);
1815  tprof_pecm_time_flatW_xview_hiW->Fill(EventTimeHigh,PE/path);
1816  }
1817  // Check non-swapped channels
1818  if( det.find("fd")!=std::string::npos && checkNoSwap ){
1819  if( checkSwap(plane,int(cell)/32) ){
1820  tprof_resp_time_flatW_xview_noSwap->Fill(EventTimeHigh,PECorr/path);
1821  tprof_pecm_time_flatW_xview_noSwap->Fill(EventTimeHigh,PE/path);
1822  }
1823  }
1824  }
1825  if(view==1 && flatW_y_min < w && w < flatW_y_max){
1826  tprof_resp_time_flatW_yview->Fill(EventTimeHigh,PECorr/path);
1827  tprof_pecm_time_flatW_yview->Fill(EventTimeHigh,PE/path);
1828  tprof_pe_time_flatW_yview->Fill(EventTimeHigh,PE);
1829  if( byDB ){
1830  tprof_pecm_time_flatW_yview_diblock[diblock-1]->Fill(EventTimeHigh,PE/path);
1831  }
1832  // W-breakdowns
1833  if( w>=flatW_y_min && w<flatW_y_min+flatW_y_selRgn ){
1834  //std::cout << "LO W = " << w << std::endl;
1835  tprof_resp_time_flatW_yview_loW->Fill(EventTimeHigh,PECorr/path);
1836  tprof_pecm_time_flatW_yview_loW->Fill(EventTimeHigh,PE/path);
1837  }
1838  else if( w>=flatW_y_min+flatW_y_selRgn && w<flatW_y_min+(2.*flatW_y_selRgn) ){
1839  //std::cout << "MD W = " << w << std::endl;
1840  tprof_resp_time_flatW_yview_mdW->Fill(EventTimeHigh,PECorr/path);
1841  tprof_pecm_time_flatW_yview_mdW->Fill(EventTimeHigh,PE/path);
1842  }
1843  else if( w>=flatW_y_min+(2.*flatW_y_selRgn) && w<flatW_y_min+(3.*flatW_y_selRgn) ){
1844  //std::cout << "HI W = " << w << std::endl;
1845  tprof_resp_time_flatW_yview_hiW->Fill(EventTimeHigh,PECorr/path);
1846  tprof_pecm_time_flatW_yview_hiW->Fill(EventTimeHigh,PE/path);
1847  }
1848  // Check non-swapped channels
1849  if( det.find("fd")!=std::string::npos && checkNoSwap ){
1850  if( checkSwap(plane,int(cell)/32) ){
1851  tprof_resp_time_flatW_yview_noSwap->Fill(EventTimeHigh,PECorr/path);
1852  tprof_pecm_time_flatW_yview_noSwap->Fill(EventTimeHigh,PE/path);
1853  }
1854  }
1855  }
1856  }
1857  //////////////////////////////////////////////////////////////////////
1858 
1859  int wbin=-1;
1860  for(unsigned int j=0; j<nWslices; j++)
1861  if( minW[j]<=w && w<=maxW[j] )
1862  wbin=j;
1863  //if(wbin<0)
1864  //std::cout << "W = " << w << " not in a bin." << std::endl;
1865 
1866  tprof_resp_fbbin->Fill(bin,PECorr/path);
1867  tprof_pecm_fbbin->Fill(bin,PE/path);
1868 
1869  tprof2_pecm_cell_plane->Fill(plane,cell,PE/path);
1870 
1871  if(view==0){
1872  if( tree_i==0 ) th1d_nxyhits_xview->Fill(nXYHits,1./double(nXYHits)); // new to look at hits/track over time
1873 
1874  th1d_nhits_resp_xview_fb[bin] ->Fill(PECorr/path);
1875  th1d_nhits_resp_xview_fb[all] ->Fill(PECorr/path);
1876  th1d_nhits_truededx_xview_fb[bin] ->Fill(trueE/path);
1877  th1d_nhits_truededx_xview_fb[all] ->Fill(trueE/path);
1878  th1d_nhits_truedetruedx_xview_fb[bin] ->Fill(trueE/true_path);
1879  th1d_nhits_truedetruedx_xview_fb[all] ->Fill(trueE/true_path);
1880  th1d_nhits_recodedx_xview_fb[bin] ->Fill(recoMeV/path);
1881  th1d_nhits_recodedx_xview_fb[all] ->Fill(recoMeV/path);
1882  th1d_nhits_pecm_xview_fb[bin] ->Fill(PE/path);
1883  th1d_nhits_pecm_xview_fb[all] ->Fill(PE/path);
1884  th1d_nhits_pe_xview_fb[bin] ->Fill(PE);
1885  th1d_nhits_pe_xview_fb[all] ->Fill(PE);
1886 
1887  tprof_resp_fbbin_xview->Fill(bin,PECorr/path);
1888  tprof_pecm_fbbin_xview->Fill(bin,PE/path);
1889  tprof_pecm_plane_xview->Fill(plane,PE/path);
1890 
1891  th1d_nhits_dcosz_xview->Fill(dcosz);
1892  th1d_nhits_dcosy_xview->Fill(dcosy);
1893  tprof_recodedx_dcosz_xview->Fill(dcosz,recoMeV/path);
1894  tprof_truededx_dcosz_xview->Fill(dcosz,trueE/path);
1895  tprof_recodedx_dcosy_xview->Fill(dcosy,recoMeV/path);
1896  tprof_truededx_dcosy_xview->Fill(dcosy,trueE/path);
1897  tprof_dx_dcosz_xview->Fill(dcosz,path);
1898  tprof_dx_dcosy_xview->Fill(dcosy,path);
1899 
1900  tprof2_resp_cell_plane_xview->Fill(plane,cell,PECorr/path);
1901  tprof2_pecm_cell_plane_xview->Fill(plane,cell,PE/path);
1902  th1d_nhits_truedx_xview->Fill(true_path);
1903  th1d_nhits_dx_xview->Fill(path);
1904  th1d_nhits_dE_xview->Fill(trueE);
1905  tprof_pecm_cell_xview ->Fill(cell, PE/path);
1906  tprof_resp_cell_xview ->Fill(cell, PECorr/path);
1907  tprof_resp_plane_xview ->Fill(plane, PECorr/path);
1908  tprof2_resp_w_cell_xview ->Fill(cell,w,PECorr/path);
1909  th2d_nhits_cell_w_xview->Fill(cell,w);
1910  th2d_nhits_dx_recodedx_xview->Fill(path,recoMeV/path);
1911  th2d_nhits_dx_truededx_xview->Fill(path,trueE/path);
1912  tprof_resp_dx_xview->Fill(path,PECorr/path);
1913  if(true_path>0.){
1914  tprof_recotruedx_truedx_xview->Fill(true_path,path/true_path);
1915  tprof_dxres_truedx_xview->Fill(true_path,(path-true_path)/true_path);
1916  tprof_dxres_w_xview->Fill(w,(path-true_path)/true_path);
1917  }
1918 
1919  th1d_nhits_w_xview_fb [bin]->Fill( w );
1920  th1d_nhits_w_xview_fb [all]->Fill( w );
1921  th2d_nhits_truededx_w_xview_fb[bin]->Fill( w, trueE/path );
1922  th2d_nhits_truededx_w_xview_fb[all]->Fill( w, trueE/path );
1923  //th2d_nhits_recodedx_w_xview_fb[bin]->Fill( w, recoMeV/path );
1924  //th2d_nhits_recodedx_w_xview_fb[all]->Fill( w, recoMeV/path );
1925  tprof_resp_w_xview_fb [bin]->Fill( w, PECorr/path );
1926  tprof_resp_w_xview_fb [all]->Fill( w, PECorr/path );
1927  if( byDB ){
1928  tprof_resp_w_xview_fb_db [bin][diblock-1]->Fill( w, PECorr/path );
1929  tprof_resp_w_xview_fb_db [all][diblock-1]->Fill( w, PECorr/path );
1930  }
1931  tprof_pecm_w_xview [bin]->Fill( w, PE/path );
1932  tprof_pecm_w_xview [all]->Fill( w, PE/path );
1933  tprof_pe_w_xview [bin]->Fill( w, PE );
1934  tprof_pe_w_xview [all]->Fill( w, PE );
1935  tprof_truededx_w_xview [bin]->Fill( w, trueE/path );
1936  tprof_truededx_w_xview [all]->Fill( w, trueE/path );
1937  tprof_recodedx_w_xview [bin]->Fill( w, recoMeV/path );
1938  tprof_recodedx_w_xview [all]->Fill( w, recoMeV/path );
1939  tprof_dx_w_xview [bin]->Fill( w, path );
1940  tprof_dx_w_xview [all]->Fill( w, path );
1941  tprof_truededx_dx_xview [bin]->Fill( path, trueE/path );
1942  tprof_truededx_dx_xview [all]->Fill( path, trueE/path );
1943  tprof_recodedx_dx_xview [bin]->Fill( path, recoMeV/path );
1944  tprof_recodedx_dx_xview [all]->Fill( path, recoMeV/path );
1945  tprof_truededx_cellpplane_xview[bin]->Fill( cells_per_plane, trueE/path );
1946  tprof_truededx_cellpplane_xview[all]->Fill( cells_per_plane, trueE/path );
1947  tprof_recodedx_cellpplane_xview[bin]->Fill( cells_per_plane, recoMeV/path );
1948  tprof_recodedx_cellpplane_xview[all]->Fill( cells_per_plane, recoMeV/path );
1949  tprof_pecorrpe_w_xview_fb [bin]->Fill( w, PECorr/PE );
1950  tprof_pecorrpe_w_xview_fb [all]->Fill( w, PECorr/PE );
1951  tprof_pecorrmev_w_xview_fb [bin]->Fill( w, PECorr/trueE );
1952  tprof_pecorrmev_w_xview_fb [all]->Fill( w, PECorr/trueE );
1953  th2d_nhits_pecm_w_xview_fb [bin]->Fill( w, PE/path );
1954  th2d_nhits_pecm_w_xview_fb [all]->Fill( w, PE/path );
1955  // NOTE NOTE NOTE NOTE NOTE NOTE NOTE
1956  // NOTE: This cut has been using w_y_cut and not w_x_cut... Keeping as is for
1957  // now and commenting out the stuff with w_x_cut for now ... BUT WHY??
1958  // NOTE NOTE NOTE NOTE NOTE NOTE NOTE
1959  if( w<w_y_cut_max && w>w_y_cut_min ){
1960  tprof_resp_pe_xview_fb [bin]->Fill( PE, PECorr/path );
1961  tprof_resp_pe_xview_fb [all]->Fill( PE, PECorr/path );
1962  tprof_pecorr_pe_xview_fb [bin]->Fill( PE, PECorr );
1963  tprof_pecorr_pe_xview_fb [all]->Fill( PE, PECorr );
1964  tprof_pecorrmev_pe_xview_fb [bin]->Fill( PE, PECorr/trueE );
1965  tprof_pecorrmev_pe_xview_fb [all]->Fill( PE, PECorr/trueE );
1966  tprof_pecorrmev_dE_xview_fb [bin]->Fill( trueE, PECorr/trueE );
1967  tprof_pecorrmev_dE_xview_fb [all]->Fill( trueE, PECorr/trueE );
1968  tprof_pecorrpe_pe_xview_fb [bin]->Fill( PE, PECorr/PE );
1969  tprof_pecorrpe_pe_xview_fb [all]->Fill( PE, PECorr/PE );
1970  tprof_truededx_pe_xview_fb [bin]->Fill( PE, trueE/path );
1971  tprof_truededx_pe_xview_fb [all]->Fill( PE, trueE/path );
1972  tprof_recodedx_pe_xview_fb [bin]->Fill( PE, recoMeV/path );
1973  tprof_recodedx_pe_xview_fb [all]->Fill( PE, recoMeV/path );
1974  tprof_truededx_pecorr_xview_fb [bin]->Fill( PECorr, trueE/path );
1975  tprof_truededx_pecorr_xview_fb [all]->Fill( PECorr, trueE/path );
1976  tprof_recodedx_pecorr_xview_fb [bin]->Fill( PECorr, recoMeV/path );
1977  tprof_recodedx_pecorr_xview_fb [all]->Fill( PECorr, recoMeV/path );
1978  tprof_pecorrpe_pecorr_xview_fb [bin]->Fill( PECorr, PECorr/PE );
1979  tprof_pecorrpe_pecorr_xview_fb [all]->Fill( PECorr, PECorr/PE );
1980  }
1981  if((cell+1)%16==0 || (cell)%16==0){
1982  tprof_resp_w_xview_edgecells_fb[bin]->Fill( w, PECorr/path );
1983  tprof_resp_w_xview_edgecells_fb[all]->Fill( w, PECorr/path );
1984  }
1985  if(wbin>=0) tprof_pe_truedE_xview_w[wbin]->Fill(w,PE);
1986  }
1987  else if (view==1) {
1988  if( tree_i==0 ) th1d_nxyhits_yview->Fill(nXYHits,1./double(nXYHits)); // new to look at hits/track over time
1989 
1990  th1d_nhits_resp_yview_fb[bin] ->Fill(PECorr/path);
1991  th1d_nhits_resp_yview_fb[all] ->Fill(PECorr/path);
1992  th1d_nhits_truededx_yview_fb[bin] ->Fill(trueE/path);
1993  th1d_nhits_truededx_yview_fb[all] ->Fill(trueE/path);
1994  th1d_nhits_truedetruedx_yview_fb[bin] ->Fill(trueE/path);
1995  th1d_nhits_truedetruedx_yview_fb[all] ->Fill(trueE/path);
1996  th1d_nhits_recodedx_yview_fb[bin] ->Fill(recoMeV/path);
1997  th1d_nhits_recodedx_yview_fb[all] ->Fill(recoMeV/path);
1998  th1d_nhits_pecm_yview_fb[bin] ->Fill(PE/path);
1999  th1d_nhits_pecm_yview_fb[all] ->Fill(PE/path);
2000  th1d_nhits_pe_yview_fb[bin] ->Fill(PE);
2001  th1d_nhits_pe_yview_fb[all] ->Fill(PE);
2002 
2003  tprof_resp_fbbin_yview->Fill(bin,PECorr/path);
2004  tprof_pecm_fbbin_yview->Fill(bin,PE/path);
2005  tprof_pecm_plane_yview->Fill(plane,PE/path);
2006 
2007  th1d_nhits_dcosz_yview->Fill(dcosz);
2008  th1d_nhits_dcosy_yview->Fill(dcosy);
2009  tprof_recodedx_dcosz_yview->Fill(dcosz,recoMeV/path);
2010  tprof_truededx_dcosz_yview->Fill(dcosz,trueE/path);
2011  tprof_recodedx_dcosy_yview->Fill(dcosy,recoMeV/path);
2012  tprof_truededx_dcosy_yview->Fill(dcosy,trueE/path);
2013  tprof_dx_dcosz_yview->Fill(dcosz,path);
2014  tprof_dx_dcosy_yview->Fill(dcosy,path);
2015 
2016  tprof2_resp_cell_plane_yview->Fill(plane,cell,PECorr/path);
2017  tprof2_pecm_cell_plane_yview->Fill(plane,cell,PE/path);
2018  th1d_nhits_truedx_yview->Fill(true_path);
2019  th1d_nhits_dx_yview->Fill(path);
2020  th1d_nhits_dE_yview->Fill(trueE);
2021  tprof_pecm_cell_yview ->Fill(cell, PE/path);
2022  tprof_resp_cell_yview ->Fill(cell, PECorr/path);
2023  tprof_resp_plane_yview ->Fill(plane, PECorr/path);
2024  tprof2_resp_cell_w_yview ->Fill(w,cell,PECorr/path);
2025  th2d_nhits_w_cell_yview ->Fill(w,cell);
2026  th2d_nhits_dx_recodedx_yview->Fill(path,recoMeV/path);
2027  th2d_nhits_dx_truededx_yview->Fill(path,trueE/path);
2028  tprof_resp_dx_yview->Fill(path,PECorr/path);
2029  if(true_path>0.){
2030  tprof_recotruedx_truedx_yview->Fill(true_path,path/true_path);
2031  tprof_dxres_truedx_yview->Fill(true_path,(path-true_path)/true_path);
2032  tprof_dxres_w_yview->Fill(w,(path-true_path)/true_path);
2033  }
2034 
2035  th1d_nhits_w_yview_fb [bin]->Fill( w );
2036  th1d_nhits_w_yview_fb [all]->Fill( w );
2037  th2d_nhits_truededx_w_yview_fb[bin]->Fill( w, trueE/path );
2038  th2d_nhits_truededx_w_yview_fb[all]->Fill( w, trueE/path );
2039  tprof_resp_w_yview_fb [bin]->Fill( w, PECorr/path );
2040  tprof_resp_w_yview_fb [all]->Fill( w, PECorr/path );
2041  if( byDB ){
2042  tprof_resp_w_yview_fb_db [bin][diblock-1]->Fill( w, PECorr/path );
2043  tprof_resp_w_yview_fb_db [all][diblock-1]->Fill( w, PECorr/path );
2044  }
2045  tprof_pecm_w_yview [bin]->Fill( w, PE/path );
2046  tprof_pecm_w_yview [all]->Fill( w, PE/path );
2047  tprof_pe_w_yview [bin]->Fill( w, PE );
2048  tprof_pe_w_yview [all]->Fill( w, PE );
2049  tprof_truededx_w_yview [bin]->Fill( w, trueE/path );
2050  tprof_truededx_w_yview [all]->Fill( w, trueE/path );
2051  tprof_recodedx_w_yview [bin]->Fill( w, recoMeV/path );
2052  tprof_recodedx_w_yview [all]->Fill( w, recoMeV/path );
2053  tprof_dx_w_yview [bin]->Fill( w, path );
2054  tprof_dx_w_yview [all]->Fill( w, path );
2055  tprof_truededx_dx_yview [bin]->Fill( path, trueE/path );
2056  tprof_truededx_dx_yview [all]->Fill( path, trueE/path );
2057  tprof_recodedx_dx_yview [bin]->Fill( path, recoMeV/path );
2058  tprof_recodedx_dx_yview [all]->Fill( path, recoMeV/path );
2059  tprof_truededx_cellpplane_yview[bin]->Fill( cells_per_plane, trueE/path );
2060  tprof_truededx_cellpplane_yview[all]->Fill( cells_per_plane, trueE/path );
2061  tprof_recodedx_cellpplane_yview[bin]->Fill( cells_per_plane, recoMeV/path );
2062  tprof_recodedx_cellpplane_yview[all]->Fill( cells_per_plane, recoMeV/path );
2063  tprof_pecorrpe_w_yview_fb [bin]->Fill( w, PECorr/PE );
2064  tprof_pecorrpe_w_yview_fb [all]->Fill( w, PECorr/PE );
2065  tprof_pecorrmev_w_yview_fb [bin]->Fill( w, PECorr/trueE );
2066  tprof_pecorrmev_w_yview_fb [all]->Fill( w, PECorr/trueE );
2067  th2d_nhits_pecm_w_yview_fb [bin]->Fill( w, PE/path );
2068  th2d_nhits_pecm_w_yview_fb [all]->Fill( w, PE/path );
2069  if( w<w_y_cut_max && w>w_y_cut_min ){
2070  tprof_resp_pe_yview_fb [bin]->Fill( PE, PECorr/path );
2071  tprof_resp_pe_yview_fb [all]->Fill( PE, PECorr/path );
2072  tprof_pecorr_pe_yview_fb [bin]->Fill( PE, PECorr );
2073  tprof_pecorr_pe_yview_fb [all]->Fill( PE, PECorr );
2074  tprof_pecorrmev_pe_yview_fb [bin]->Fill( PE, PECorr/trueE );
2075  tprof_pecorrmev_pe_yview_fb [all]->Fill( PE, PECorr/trueE );
2076  tprof_pecorrmev_dE_yview_fb [bin]->Fill( trueE, PECorr/trueE );
2077  tprof_pecorrmev_dE_yview_fb [all]->Fill( trueE, PECorr/trueE );
2078  tprof_pecorrpe_pe_yview_fb [bin]->Fill( PE, PECorr/PE );
2079  tprof_pecorrpe_pe_yview_fb [all]->Fill( PE, PECorr/PE );
2080  tprof_truededx_pe_yview_fb [bin]->Fill( PE, trueE/path );
2081  tprof_truededx_pe_yview_fb [all]->Fill( PE, trueE/path );
2082  tprof_recodedx_pe_yview_fb [bin]->Fill( PE, recoMeV/path );
2083  tprof_recodedx_pe_yview_fb [all]->Fill( PE, recoMeV/path );
2084  tprof_truededx_pecorr_yview_fb [bin]->Fill( PECorr, trueE/path );
2085  tprof_truededx_pecorr_yview_fb [all]->Fill( PECorr, trueE/path );
2086  tprof_recodedx_pecorr_yview_fb [bin]->Fill( PECorr, recoMeV/path );
2087  tprof_recodedx_pecorr_yview_fb [all]->Fill( PECorr, recoMeV/path );
2088  tprof_pecorrpe_pecorr_yview_fb [bin]->Fill( PECorr, PECorr/PE );
2089  tprof_pecorrpe_pecorr_yview_fb [all]->Fill( PECorr, PECorr/PE );
2090  }
2091  if((cell+1)%16==0 || (cell)%16==0){
2092  tprof_resp_w_yview_edgecells_fb[bin]->Fill( w, PECorr/path );
2093  tprof_resp_w_yview_edgecells_fb[all]->Fill( w, PECorr/path );
2094  }
2095 
2096  if(wbin>=0) tprof_pe_truedE_yview_w[wbin]->Fill(w,PE);
2097  }
2098  else { std::cout << "unknown view" << std::endl; }
2099  }
2100  } // ttree loop
2101 
2102  // NEW: Close the input file and cd to the output file
2103  // NOTE: Closing didn't work... Close it after in the .cxx
2104  //f->Close();
2105  outFile->cd();
2106 
2107  std::cout << "# hits passing track window, w, and cell cuts: " << n_track_window_hits << std::endl;
2108  std::cout << "# hits which were not from muons: " << n_nonmuon_hits << std::endl;
2109 
2110  // Print some output on this run
2111  if(n_zero_true_path!=0)
2112  std::cout << "truepath was nonpositive " << n_zero_true_path << " times" << std::endl;
2113  std::cout << "path was nonpositive " << n_zero_path << " times" << std::endl;
2114  std::cout << "PECorr was zero " << n_zero_PEc << " times" << std::endl;
2115  std::cout << "PECorr was negative " << n_neg_PEc << " times" << std::endl;
2116  std::cout << "trueE was nonpositive " << n_zero_trueE << " times" << std::endl;
2117 
2118  metadataPlot->Write();
2119 
2120  th1d_nhits_resp_flatW->Write();
2121  th1d_nhits_resp_flatW_xview->Write();
2122  th1d_nhits_resp_flatW_yview->Write();
2123  th1d_nhits_recodedx_flatW_xview->Write();
2124  th1d_nhits_recodedx_flatW_yview->Write();
2125 
2126  ///// Additional plots
2127  if( nDivisions>0 ){
2128  for( unsigned int thisTime=0; thisTime<nDivisions; ++thisTime ){
2129  th1d_nhits_resp_flatW_range[thisTime]->Write();
2130  th1d_nhits_resp_flatW_xview_range[thisTime]->Write();
2131  th1d_nhits_resp_flatW_yview_range[thisTime]->Write();
2132  }
2133  }
2134  if( byDB ){
2135  for( unsigned int iChan=0; iChan<nChan; ++iChan ){
2136  th1d_nhits_resp_flatW_xview_diblock[iChan]->Write();
2137  th1d_nhits_resp_flatW_xview_diblock[iChan]->Write();
2138  th1d_nhits_resp_flatW_yview_diblock[iChan]->Write();
2139  if( nDivisions>0 ){
2140  for( unsigned int thisTime=0; thisTime<nDivisions; ++thisTime ){
2141  th1d_nhits_resp_flatW_range_diblock[thisTime][iChan]->Write();
2142  th1d_nhits_resp_flatW_xview_range_diblock[thisTime][iChan]->Write();
2143  th1d_nhits_resp_flatW_yview_range_diblock[thisTime][iChan]->Write();
2144  }
2145  }
2146  }
2147  }
2148  //////////////////////
2149 
2150  tprof_pecm_plane_xview->Write();
2151  tprof_pecm_plane_yview->Write();
2152 
2153  th1d_nhits_dcosz_xview->Write(); th1d_nhits_dcosz_yview->Write();
2154  th1d_nhits_dcosy_xview->Write(); th1d_nhits_dcosy_yview->Write();
2155  tprof_recodedx_dcosz_xview->Write(); tprof_recodedx_dcosz_yview->Write();
2156  tprof_truededx_dcosz_xview->Write(); tprof_truededx_dcosz_yview->Write();
2157  tprof_recodedx_dcosy_xview->Write(); tprof_recodedx_dcosy_yview->Write();
2158  tprof_truededx_dcosy_xview->Write(); tprof_truededx_dcosy_yview->Write();
2159  tprof_dx_dcosz_xview->Write(); tprof_dx_dcosz_yview->Write();
2160  tprof_dx_dcosy_xview->Write(); tprof_dx_dcosy_yview->Write();
2161 
2162  th1d_nhits_dx_xview->Write();
2163  th1d_nhits_dx_yview->Write();
2164  tprof2_resp_cell_plane_xview->Write();
2165  tprof2_resp_cell_plane_yview->Write();
2166  tprof_pecm_cell_xview->Write();
2167  tprof_pecm_cell_yview->Write();
2168  tprof_resp_cell_xview->Write();
2169  tprof_resp_cell_yview->Write();
2170  tprof_resp_plane_xview->Write();
2171  tprof_resp_plane_yview->Write();
2172  tprof2_resp_w_cell_xview->Write();
2173  tprof2_resp_cell_w_yview->Write();
2174  th2d_nhits_cell_w_xview->Write();
2175  th2d_nhits_w_cell_yview->Write();
2176  th2d_nhits_dx_recodedx_xview->Write();
2177  th2d_nhits_dx_truededx_xview->Write();
2178  th2d_nhits_dx_recodedx_yview->Write();
2179  th2d_nhits_dx_truededx_yview->Write();
2180  tprof2_pecm_cell_plane_xview->Write();
2181  tprof2_pecm_cell_plane_yview->Write();
2182  tprof2_pecm_cell_plane->Write();
2183  th2d_uncalib_cell_plane->Write(); // uncalibrated cells plot
2184  th2d_uncalib_cell_plane_xview->Write();
2185  th2d_uncalib_cell_plane_yview->Write();
2186  tprof_resp_fbbin->Write();
2187  tprof_pecm_fbbin->Write();
2188  tprof_resp_fbbin_xview->Write();
2189  tprof_pecm_fbbin_xview->Write();
2190  tprof_resp_fbbin_yview->Write();
2191  tprof_pecm_fbbin_yview->Write();
2192  tprof_resp_dx_xview->Write();
2193  tprof_resp_dx_yview->Write();
2194 
2195  if(isMC){
2196  th2d_nhits_dE_pe_xview->Write();
2197  th2d_nhits_dE_pe_yview->Write();
2198  th1d_nhits_truedx_xview->Write();
2199  th1d_nhits_truedx_yview->Write();
2200  th1d_nhits_truededx_flatW->Write();
2201  th1d_nhits_truededx_flatW_xview->Write();
2202  th1d_nhits_truededx_flatW_yview->Write();
2203  ///// Additional plots
2204  if( byDB ){
2205  for( unsigned int iChan=0; iChan<nChan; ++iChan ){
2206  th1d_nhits_truededx_flatW_diblock[iChan]->Write();
2207  th1d_nhits_truededx_flatW_xview_diblock[iChan]->Write();
2208  th1d_nhits_truededx_flatW_yview_diblock[iChan]->Write();
2209  }
2210  }
2211  //////////////////////
2212  tprof_recotruedx_truedx_xview->Write();
2213  tprof_recotruedx_truedx_yview->Write();
2214  tprof_dxres_truedx_xview->Write();
2215  tprof_dxres_truedx_yview->Write();
2216  tprof_dxres_w_xview->Write();
2217  tprof_dxres_w_yview->Write();
2218  th1d_nhits_dE_xview->Write();
2219  th1d_nhits_dE_yview->Write();
2220  }
2221 
2222  for(unsigned int b=0; b<nFBbins+1; b++){
2223  // vs W
2224  th1d_nhits_w_xview_fb[b]->Write();
2225  th1d_nhits_w_yview_fb[b]->Write();
2226  tprof_resp_w_xview_edgecells_fb[b]->Write();
2227  tprof_resp_w_yview_edgecells_fb[b]->Write();
2228  tprof_resp_w_xview_fb[b]->Write();
2229  tprof_resp_w_yview_fb[b]->Write();
2230  if( byDB ){
2231  for( unsigned int thisDB=0; thisDB<nChan; ++thisDB ){
2232  tprof_resp_w_xview_fb_db[b][thisDB]->Write();
2233  tprof_resp_w_yview_fb_db[b][thisDB]->Write();
2234  }
2235  }
2236  tprof_pecm_w_xview[b]->Write();
2237  tprof_pecm_w_yview[b]->Write();
2238  tprof_recodedx_w_xview[b]->Write();
2239  tprof_recodedx_w_yview[b]->Write();
2240  tprof_pe_w_xview[b]->Write();
2241  tprof_pe_w_yview[b]->Write();
2242  tprof_dx_w_xview[b]->Write();
2243  tprof_dx_w_yview[b]->Write();
2244  tprof_pecorrpe_w_xview_fb[b]->Write();
2245  tprof_pecorrpe_w_yview_fb[b]->Write();
2246  // vs PE
2247  tprof_resp_pe_xview_fb[b]->Write();
2248  tprof_resp_pe_yview_fb[b]->Write();
2249  tprof_pecorr_pe_xview_fb[b]->Write();
2250  tprof_pecorr_pe_yview_fb[b]->Write();
2251  tprof_pecorrpe_pe_xview_fb[b]->Write();
2252  tprof_pecorrpe_pe_yview_fb[b]->Write();
2253  tprof_recodedx_pe_xview_fb[b]->Write();
2254  tprof_recodedx_pe_yview_fb[b]->Write();
2255  // vs PECor
2256  tprof_recodedx_pecorr_xview_fb[b]->Write();
2257  tprof_recodedx_pecorr_yview_fb[b]->Write();
2258  tprof_pecorrpe_pecorr_xview_fb[b]->Write();
2259  tprof_pecorrpe_pecorr_yview_fb[b]->Write();
2260  // vs x
2261  tprof_recodedx_x_xview[b]->Write();
2262  tprof_recodedx_x_yview[b]->Write();
2263  // vs dx, pathlength
2264  tprof_recodedx_dx_xview[b]->Write();
2265  tprof_recodedx_dx_yview[b]->Write();
2266  // vs cells per plane
2267  tprof_recodedx_cellpplane_xview[b]->Write();
2268  tprof_recodedx_cellpplane_yview[b]->Write();
2269  // vs FB
2270  th1d_nhits_resp_xview_fb[b]->Write();
2271  th1d_nhits_resp_yview_fb[b]->Write();
2272  th1d_nhits_recodedx_xview_fb[b]->Write();
2273  th1d_nhits_recodedx_yview_fb[b]->Write();
2274  th1d_nhits_pecm_xview_fb[b]->Write();
2275  th1d_nhits_pecm_yview_fb[b]->Write();
2276  th1d_nhits_pe_xview_fb[b]->Write();
2277  th1d_nhits_pe_yview_fb[b]->Write();
2278  // Other
2279  th2d_nhits_pecm_w_xview_fb[b]->Write();
2280  th2d_nhits_pecm_w_yview_fb[b]->Write();
2281 
2282  if(isMC){
2283  th2d_nhits_truededx_w_xview_fb[b]->Write();
2284  th2d_nhits_truededx_w_yview_fb[b]->Write();
2285  th2d_nhits_truededx_w_xview_nowindow_fb[b]->Write();
2286  th2d_nhits_truededx_w_yview_nowindow_fb[b]->Write();
2287  tprof_truededx_w_xview[b]->Write();
2288  tprof_truededx_w_yview[b]->Write();
2289  tprof_pecorrmev_w_xview_fb[b]->Write();
2290  tprof_pecorrmev_w_yview_fb[b]->Write();
2291  tprof_pecorrmev_pe_xview_fb[b]->Write();
2292  tprof_pecorrmev_pe_yview_fb[b]->Write();
2293  tprof_truededx_pe_xview_fb[b]->Write();
2294  tprof_truededx_pe_yview_fb[b]->Write();
2295  tprof_truededx_pecorr_xview_fb[b]->Write();
2296  tprof_truededx_pecorr_yview_fb[b]->Write();
2297  tprof_pecorrmev_dE_xview_fb[b]->Write();
2298  tprof_pecorrmev_dE_yview_fb[b]->Write();
2299  tprof_truededx_x_xview[b]->Write();
2300  tprof_truededx_x_yview[b]->Write();
2301  tprof_truededx_dx_xview[b]->Write();
2302  tprof_truededx_dx_yview[b]->Write();
2303  tprof_truededx_cellpplane_xview[b]->Write();
2304  tprof_truededx_cellpplane_yview[b]->Write();
2305  th1d_nhits_truededx_xview_fb[b]->Write();
2306  th1d_nhits_truededx_yview_fb[b]->Write();
2307  th1d_nhits_truedetruedx_xview_fb[b]->Write();
2308  th1d_nhits_truedetruedx_yview_fb[b]->Write();
2309  }
2310  }
2311 
2312  for(unsigned int wb=0; wb<nWslices; wb++){
2313  if(isMC){
2314  tprof_pe_truedE_xview_w[wb]->Write();
2315  tprof_pe_truedE_yview_w[wb]->Write();
2316  }
2317  }
2318 
2319  if( tree_i==0 ){
2320  tprof_resp_time_flatW_xview->Write();
2321  tprof_pecm_time_flatW_xview->Write();
2322  tprof_pe_time_flatW_xview->Write();
2323  tprof_resp_time_flatW_yview->Write();
2324  tprof_pecm_time_flatW_yview->Write();
2325  tprof_pe_time_flatW_yview->Write();
2326 
2327  tprof_resp_time_flatW_xview_loW->Write();
2328  tprof_pecm_time_flatW_xview_loW->Write();
2329  tprof_resp_time_flatW_yview_loW->Write();
2330  tprof_pecm_time_flatW_yview_loW->Write();
2331  tprof_resp_time_flatW_xview_mdW->Write();
2332  tprof_pecm_time_flatW_xview_mdW->Write();
2333  tprof_resp_time_flatW_yview_mdW->Write();
2334  tprof_pecm_time_flatW_yview_mdW->Write();
2335  tprof_resp_time_flatW_xview_hiW->Write();
2336  tprof_pecm_time_flatW_xview_hiW->Write();
2337  tprof_resp_time_flatW_yview_hiW->Write();
2338  tprof_pecm_time_flatW_yview_hiW->Write();
2339 
2340  if( det.find("fd")!=std::string::npos && checkNoSwap ){
2341  tprof_resp_time_flatW_xview_noSwap->Write();
2342  tprof_pecm_time_flatW_xview_noSwap->Write();
2343  tprof_resp_time_flatW_yview_noSwap->Write();
2344  tprof_pecm_time_flatW_yview_noSwap->Write();
2345  }
2346 
2347  if( byDB ){
2348  for(unsigned int iChan=0; iChan<tprof_pecm_time_flatW_xview_diblock.size(); ++iChan){
2349  tprof_pecm_time_flatW_xview_diblock[iChan]->Write();
2350  tprof_pecm_time_flatW_yview_diblock[iChan]->Write();
2351  }
2352  }
2353 
2354  th1d_nxyhits_xview->Write();
2355  th1d_nxyhits_yview->Write();
2356  th1d_nxyhits_flatW_xview->Write();
2357  th1d_nxyhits_flatW_yview->Write();
2358  }
2359 
2360  outFile->Close();
2361 
2362  return;
2363 }
const XML_Char * name
Definition: expat.h:151
diblock
print "ROW IS " print row
Definition: geo2elec.py:31
void RunCalibration(TFile *f, std::string det, bool isMC, std::string inFileName, std::string outFileName, long maxEntries, bool byDB, int runsPerDiv, bool checkNoSwap)
Definition: Fennec.h:307
T sqrt(T number)
Definition: d0nt_math.hpp:156
constexpr T pow(T x)
Definition: pow.h:75
cout<< t1-> GetEntries()
Definition: plottest35.C:29
float abs(float number)
Definition: d0nt_math.hpp:39
TFile * outFile
Definition: PlotXSec.C:135
bool checkSwap(int plane, int mod)
Definition: FennecFnc.h:83
correl_xv GetXaxis() -> SetDecimals()
inFileName
if we get here, we&#39;re doing the base definitions:
Definition: mkDefs.py:168
const double j
Definition: BetheBloch.cxx:29
float bin[41]
Definition: plottest35.C:14
OStream cout
Definition: OStream.cxx:6
const std::string path
Definition: plot_BEN.C:43
std::string getOutputName(std::string detector, std::string period, std::string sample, std::string outName, bool appendRoot, std::string user_dir)
Function to determine output file name.
Definition: FennecFnc.h:221
const hit & b
Definition: hits.cxx:21
void MakeCSV(TFile *inFileData, TFile *inFileMC, std::string detector, std::string sample, std::string period, std::string outputName, bool byDB, int runsPerDiv, std::string calibMethod, bool mcRun, std::string user_dir, bool verbosity)
Definition: Fennec.h:4
void Format(TGraph *gr, int lcol, int lsty, int lwid, int mcol, int msty, double msiz)
Definition: Style.cxx:154
std::string to_string(ModuleType mt)
Definition: ModuleType.h:32
Float_t e
Definition: plot.C:35
Float_t w
Definition: plot.C:20
bool cut_Track_Window(double x)
Definition: FennecFnc.h:100
fvar< T > ceil(const fvar< T > &x)
Definition: ceil.hpp:11
bool cut_PE(double pe)
Definition: FennecFnc.h:104
enum BeamMode string