20 #include "TFitResult.h"
26 using namespace KM3NETDAQ;
28 int main(
int argc ,
char** argv){
34 bool overwriteDetector;
45 catch(
const exception &error) {
46 ERROR(error.what() << endl);
50 load(detectorFile , detector);
57 TFile
in(inputFiles.c_str() ,
"read");
61 string outFile_root = outFile+
".root";
62 string outFile_csv = outFile+
".csv";
63 string outFile_det = outFile+
"_det.csv";
65 TFile output(outFile_root.c_str(),
"RECREATE");
66 cout <<
" Writing histograms to: " << outFile_root << endl;
69 ofstream out_csv(outFile_csv.c_str());
70 cout <<
" Writing fitting data to: " << outFile_csv << endl;
71 out_csv <<
"String,Ref_floor,Tgt_floor,T0,Error" << endl;
79 for(
int j = 1;
j<=num_of_strings;
j++){
80 vmap.insert( make_pair(
j,base_vector) );
83 for(JDetector::iterator module = detector.begin(); module != detector.end(); ++module){
84 int string_no = moduleRouter.
getModule( module->getID() ).getString();
85 int tgt_floor = moduleRouter.
getModule( module->getID() ).getFloor();
90 for(
int i = 1; i <= Neighbours; i++){
91 int ref_floor = tgt_floor - i;
93 if(ref_floor < 1){
continue; }
97 if(
in.GetListOfKeys()->Contains(hist_str.c_str()) ){
99 in.GetObject( hist_str.c_str(), h2D );
102 int bin = c.
getIndex(tgt_floor - 1, ref_floor - 1 ) + 1;
104 hist_name =
"S"+hist_str+
"F"+
to_string(ref_floor)+
"_S"+hist_str+
"F"+
to_string(tgt_floor);
105 TH1D* h1D = h2D->ProjectionY( hist_name.c_str(), bin, bin );
106 if( h1D->GetEntries() < 1 ){
continue; }
109 T0s.push_back( fitresult->Parameter(1) );
110 Errors.push_back( fitresult->ParError(1) );
113 out_csv << string_no <<
"," << ref_floor <<
"," << tgt_floor <<
"," << fitresult->Parameter(1) <<
"," << fitresult->ParError(1) << endl;
115 if(!outFile.empty()){
122 double weighted_T0 = 0;
124 double Sum_of_Weight = 0;
125 for(
size_t i = 0, max = T0s.size(); i != max; ++i ){
126 Weight = 1/Errors[i];
127 weighted_T0 += T0s[i]*Weight;
128 Sum_of_Weight += Weight;
130 if(Sum_of_Weight != 0){
131 weighted_T0 /= Sum_of_Weight;
135 map_type::iterator map_iterator = vmap.find(string_no);
136 map_iterator->second[tgt_floor-1] = weighted_T0;
141 for(
int j = 1;
j<=num_of_strings;
j++){
142 map_type::iterator map_iterator = vmap.find(
j);
144 int n = map_iterator->second.size();
145 double average =
accumulate( map_iterator->second.begin(), map_iterator->second.end(), 0.0)/
n;
147 for(
int j = 0;
j<
n;
j++ ){
148 map_iterator->second[
j] = map_iterator->second[
j] - average;
152 ofstream out_det(outFile_det.c_str());
153 cout <<
" Writing det changes to: " << outFile_det << endl;
154 out_det <<
"String,Floor,T0" << endl;
157 for(JDetector::iterator module = detector.begin(); module != detector.end(); ++module){
158 int string_no = moduleRouter.
getModule( module->getID() ).getString();
159 int tgt_floor = moduleRouter.
getModule( module->getID() ).getFloor();
161 map_type::iterator map_iterator = vmap.find(string_no);
162 double weighted_T0 = map_iterator->second[tgt_floor-1];
164 if(overwriteDetector){
169 out_det << string_no <<
"," << tgt_floor <<
"," << weighted_T0 << endl;
173 if(overwriteDetector){
174 cout <<
" Overwriting the detectorfile" << endl;
176 store(detectorFile,detector);
Utility class to parse command line options.
const JModule & getModule(const JObjectID &id) const
Get module parameters.
TFitResultPtr Fit(TH1D *h)
Router for direct addressing of module data in detector data structure.
Data structure for detector geometry and calibration.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
void store(const std::string &file_name, const JDetector &detector)
Store detector to output file.
int getNumberOfFloors(const JDetector &detector)
Get number of floors.
const JPMT & getPMT(const int index) const
Get PMT.
static const JStringCounter getNumberOfStrings
Function object to count unique strings.
Direct access to module in detector data structure.
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Utility class to parse command line options.
alias put_queue eval echo n
std::string to_string(const T &value)
Convert value to string.
void addT0(const double t0)
Add time offset.
bool comparepair(const pair_type &A, const pair_type &B)
do set_variable DETECTOR_TXT $WORKDIR detector
KM3NeT DAQ constants, bit handling, etc.
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY source JAcoustics sh $DETECTOR_ID typeset A TRIPODS get_tripods $WORKDIR tripod txt TRIPODS for EMITTER in
int main(int argc, char *argv[])