Jpp
Classes | Functions
JVoltageOptimizer_IO.hh File Reference
#include <iostream>
#include "JDAQ/JDAQTimeslice.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JModuleLocation.hh"

Go to the source code of this file.

Classes

struct  IO
 Structure to store the different command line arguments for JRunAnalyzer. More...
 

Functions

int read_user_options (IO &options, int argc, char **argv)
 Parses the command line options and fills an IO structure with them. More...
 
void write_output_opt (string filename, TTree *t, vector< TH2D * > summary_refs, vector< NBRun * > &Runs)
 Writes a .root file with the optional output produced by JVoltageOptimizer. More...
 
void write_output_txt (string filename, vector< pair< double, int > > optimal_voltages)
 Writes a .txt file with the result produced by JVoltageOptimizer. More...
 
void write_output_cal (string filename, vector< pair< double, int > > &optimal_voltages, vector< NBRun * > &Runs, int string_number, int nFloors)
 Writes a .root file with the optional output produced by JVoltageOptimizer. More...
 

Function Documentation

◆ read_user_options()

int read_user_options ( IO options,
int  argc,
char **  argv 
)
inline

Parses the command line options and fills an IO structure with them.

Parameters
optionsan option structure
argcthe number of command line arguments
argvthe command line arguments
Returns
1 if works 2 if it doesn't work

Definition at line 68 of file JVoltageOptimizer_IO.hh.

68  {
69 
70  int a = 0;
71 
72  try {
73 
74  JParser<string> zap;
75 
76  zap["o"] = make_field(options.ofname_opt) = "out_opt.root" ;
77 
78  zap["c"] = make_field(options.ofname_cal) = "out_cal.root" ;
79 
80  zap["t"] = make_field(options.ofname_txt) = "out.txt" ;
81 
82  zap["f"] = make_field(options.ifnames) ;
83 
84  zap["a"] = make_field(options.detector_file) ;
85 
86  zap["s"] = make_field(options.string_number) ;
87 
88  zap["u"] = make_field(options.up_pmts) = 1 ;
89 
90  zap["d"] = make_field(options.down_pmts) = 4 ;
91 
92  zap["n"] = make_field(options.max_neighbors) = 2 ;
93 
94  zap["S"] = make_field(options.saturation_tot) = 26.4 ;
95 
96  if (zap.read(argc, argv) != 0)
97 
98  a = 1;
99 
100  }
101 
102  catch(const exception &error) {
103 
104  ERROR(error.what() << endl);
105 
106  a = 2;
107 
108  }
109 
110  return a;
111 
112 }

◆ write_output_opt()

void write_output_opt ( string  filename,
TTree *  t,
vector< TH2D * >  summary_refs,
vector< NBRun * > &  Runs 
)
inline

Writes a .root file with the optional output produced by JVoltageOptimizer.

Parameters
filenameto output .root filex
ta root TTree.
summary_refsa TH2 with information of which reference pmts were used for optimizing the nano-beacon voltages.
RunsThe list of nanobeacon Runs analyzed for the optimization.

Definition at line 123 of file JVoltageOptimizer_IO.hh.

123  {
124 
125  TFile* f = new TFile(filename.c_str() , "recreate");
126 
127  f->cd();
128 
129  t->Write();
130 
131  int i = 0 ;
132 
133  char name [100] ;
134 
135  for (auto & run : Runs){
136 
137  sprintf(name , "RUN_%d_%2.1fV" , run->getRunNumber() , run->getVoltage()) ;
138 
139  f->mkdir(name) ;
140 
141  f->cd(name) ;
142 
143  summary_refs[i]->Write() ;
144 
145  for(auto & sm : run->getSuperModules()){
146 
147  sprintf(name , "RUN_%d_%2.1fV/NB_%d" , run->getRunNumber() , run->getVoltage() , sm->getFloor()) ;
148 
149  f->mkdir(name) ;
150 
151  f->cd(name) ;
152 
153  for(auto & spm : sm->get_ref_pmts()){
154 
155  spm->getNBPulse()->gettime_vs_tot()->Write() ;
156 
157  spm->getNBPulse()->gettime_peak()->Write() ;
158 
159  spm->getNBPulse()->gettot_peak()->Write() ;
160 
161  }
162 
163  }
164 
165  i++ ;
166 
167  }
168 
169  f->Close() ;
170 
171  delete(f) ;
172 
173 }

◆ write_output_txt()

void write_output_txt ( string  filename,
vector< pair< double, int > >  optimal_voltages 
)
inline

Writes a .txt file with the result produced by JVoltageOptimizer.

This file contains a list of nanobeacons and optimal voltages

Parameters
filenamepath to the name of the .txt file
optimal_voltagesList of pairs containing optimal voltage and the corresponding run number for each nanobeacon

Definition at line 183 of file JVoltageOptimizer_IO.hh.

183  {
184 
185  ofstream txtfile;
186 
187  txtfile.open (filename);
188 
189  txtfile << "# Table of optimal voltages: column 1 = floor , column 2 = voltage" << endl;
190 
191  int i=0 ;
192 
193  for (auto & op : optimal_voltages){
194 
195  txtfile << i+1 << "\t" << op.first << endl;
196 
197  i++ ;
198 
199  }
200 
201  txtfile.close();
202 
203 }

◆ write_output_cal()

void write_output_cal ( string  filename,
vector< pair< double, int > > &  optimal_voltages,
vector< NBRun * > &  Runs,
int  string_number,
int  nFloors 
)
inline

Writes a .root file with the optional output produced by JVoltageOptimizer.

This file has the format needed by JInterDomCal to perform a calibration of the DU and contains the pulses corresponding to the optimized voltages.

Parameters
filenamepath to the name of the .root
optimal_voltagesList of pairs containing optimal voltage and the corresponding run number for each nanobeacon
Runsvector of NBruns
string_numberNumber of the DU to calibrate.
nFloorsNumber of floors in the detector.

Definition at line 215 of file JVoltageOptimizer_IO.hh.

215  {
216 
217  TFile* out_file = new TFile(filename.c_str() , "recreate") ;
218 
219  TVectorD Run (2) ;
220 
221  Run[0] = 999 ;
222 
223  Run[1] = 999 ;
224 
225  Run.Write("Run_Info") ;
226 
227  int Ndoms = nFloors ;
228 
229  char dirname[200] ;
230 
231  for( int i = 0 ; i < Ndoms ; ++i ) {
232 
233  sprintf(dirname, "S%iF%i", string_number, i + 1 ) ;
234 
235  out_file->mkdir(dirname)->cd() ;
236 
237  }
238 
239  out_file->cd();
240 
241  int nb_index = 0;
242 
243  for (auto & volt :optimal_voltages){
244 
245  NBRun* r = Runs[volt.second] ;
246 
247  vector< pair < SuperModule* , vector < SuperPMT* > > > targets = r->getSuperModules()[nb_index]->get_targets() ;
248 
249  for (auto & tgt : targets){
250 
251  int floor = tgt.first->getFloor() ;
252 
253  sprintf( dirname , "S%iF%i", string_number , floor ) ;
254 
255  out_file->cd(dirname) ;
256 
257  for (auto & spm : tgt.second){
258 
259  spm->getNBPulse()->gettime_vs_tot()->Write() ;
260 
261  }
262 
263  }
264 
265  sprintf(dirname , "S%iF%i", string_number , r->getSuperModules()[nb_index]->getFloor()) ;
266 
267  out_file->cd(dirname) ;
268 
269  for (auto & spm : r->getSuperModules()[nb_index]->get_ref_pmts()){
270 
271  spm->getNBPulse()->gettime_vs_tot()->Write() ;
272 
273  }
274 
275  out_file->cd();
276 
277  nb_index ++ ;
278 
279  }
280 
281  out_file->Close();
282 
283  delete(out_file) ;
284 
285 }
IO::max_neighbors
int max_neighbors
Definition: JVoltageOptimizer_IO.hh:51
IO::ofname_opt
string ofname_opt
Definition: JVoltageOptimizer_IO.hh:37
IO::ofname_cal
string ofname_cal
Definition: JVoltageOptimizer_IO.hh:39
IO::ifnames
JMultipleFileScanner ifnames
Definition: JPulseFinder_IO.hh:40
IO::down_pmts
int down_pmts
Definition: JInterDomCal_IO.hh:59
NBRun
Class dedicated to the nanobeacon analyses, where the Modules in the detector are not regarded as sin...
Definition: NBRun.hh:24
std::vector
Definition: JSTDTypes.hh:12
IO::saturation_tot
double saturation_tot
Definition: JVoltageOptimizer_IO.hh:53
JPARSER::JParser
Utility class to parse command line options.
Definition: JParser.hh:1493
ERROR
#define ERROR(A)
Definition: JMessage.hh:66
IO::up_pmts
int up_pmts
Definition: JInterDomCal_IO.hh:57
IO::ofname_txt
string ofname_txt
Definition: JInterDomCal_IO.hh:51
IO::string_number
int string_number
Definition: JInterDomCal_IO.hh:55
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1954
IO::detector_file
string detector_file
Definition: JRunAnalyzer.cc:46
JTOOLS::r
data_type r[M+1]
Definition: JPolint.hh:709
JPARSER::JParser::read
int read(const int argc, const char *const argv[])
Parse the program's command line options.
Definition: JParser.hh:1791