Jpp 19.3.0-rc.2
the software that should make you happy
Loading...
Searching...
No Matches
JOffset_extract.cc
Go to the documentation of this file.
1\
2#include <iostream>
3#include <string>
4#include <map>
5
6#include "TFile.h"
7#include "TH1D.h"
8#include "TH2D.h"
9#include "TF1.h"
10#include "TMath.h"
11
12#include "JROOT/JMinimizer.hh"
13
17
18#include "Jeep/JPrint.hh"
19#include "Jeep/JParser.hh"
20
21
22/**
23 * \file
24 *
25 * Auxiliary program to extract string-string correlation information from the output created with JMonitorL1dt
26 *
27 * \author dsamtleben
28 */
29int main(int argc, char **argv)
30{
31 using namespace std;
32 using namespace JPP;
33
34 string inputFile;
35 string outputFile;
36 string detectorFile;
37 int neighbour, setmax;
38 int debug;
39
40 try {
41
42 JParser<> zap("Program to extract time offsets of DOM-DOM correlations");
43
44 zap['f'] = make_field(inputFile, "input file") = "monitor.root";
45 zap['o'] = make_field(outputFile, "output file") = "corr_max.root";
46 zap['a'] = make_field(detectorFile, "detector file");
47 zap['n'] = make_field(neighbour, "neighbour level") = 5;
48 zap['m'] = make_field(setmax, "minimal entries") = 5;
49 zap['d'] = make_field(debug) = 0;
50
51 if (zap.read(argc, argv) != 0) {
52 return 1;
53 }
54 }
55 catch(const exception &error) {
56 FATAL(error.what() << endl);
57 }
58
59
61
62 try {
63 load(detectorFile, detector);
64 }
65 catch(const JException& error) {
66 FATAL(error);
67 }
68
69
70 const JStringRouter string(detector);
71
73
74 for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
75 zmap[module->getString()][module->getFloor()] = module->getID();
76 }
77
78 struct h1_t {
79 h1_t() :
80 p(NULL) // pointer must be initialised
81 {}
82
83 TH1D* p;
84 };
85
86 TF1 f1("f1", "[0]*exp(-0.5*(x-[1])*(x-[1])/([2]*[2])) + [3]");
87
88 const int number_of_strings = getNumberOfStrings(detector);
89
90// cout << "number of strings " << number_of_strings << endl;
91
92 double maxarr[number_of_strings][number_of_strings][3];
93 int duarr[number_of_strings];
94
96
97 TFile* in = TFile::Open(inputFile.c_str(), "exist");
98
99 if (in == NULL || !in->IsOpen()) {
100 FATAL("File: " << inputFile << " not opened." << endl);
101 }
102
103
104 for (const auto& string_1 : zmap) {
105 duarr[string.getIndex(string_1.first)]=string_1.first;
106
107 for (const auto& floor_1 : string_1.second) {
108
109 const int module_1 = floor_1.second; // parent module
110
111 TH2D* h2 = (TH2D*) in->Get(MAKE_CSTRING(module_1 << ".2S"));
112
113 if (h2 == NULL) {
114 continue;
115 }
116
117 for (const auto& string_2 : zmap) {
118 duarr[string.getIndex(string_2.first)]=string_2.first;
119
120 if (string_1.first != string_2.first) {
121
122 TH1D* h1 = H1[string_1.first][string_2.first].p;
123
124 if (h1 == NULL) {
125
126 h1 = new TH1D(MAKE_CSTRING(string_1.first << "_" << string_2.first << "_" << neighbour << ".2T"), NULL,
127 h2->GetYaxis()->GetNbins(), h2->GetYaxis()->GetXmin(), h2->GetYaxis()->GetXmax());
128
129 H1[string_1.first][string_2.first].p = h1; // book keeping
130 }
131
132 for (const auto& floor_2 : string_2.second) {
133 if (floor_1.first > floor_2.first && (floor_1.first - floor_2.first) == neighbour) { // floor condition
134
135 const int module_2 = floor_2.second; // daughter module
136
137 TH1D* py = h2->ProjectionY("__py",
138 h2->GetXaxis()->FindBin(TString(Form("%i", module_2))),
139 h2->GetXaxis()->FindBin(TString(Form("%i", module_2))), "e");
140
141 h1->Add(py); // add data
142
143 delete py;
144 }
145 }
146//------------
147 if (floor_1.first == 18){
148
149 double mm=h1->GetXaxis()->GetBinCenter(h1->GetMaximumBin());
150
151 f1.SetParameter(0, double(h1->GetMaximum()));
152 f1.SetParameter(1, double(mm));
153 f1.SetParameter(2, 50.);
154
155 for (Int_t i = 0; i != f1.GetNpar(); ++i) {
156 f1.SetParError(i, 0.0);
157 }
158
159 // fit
160
161 if (h1->GetMaximum()>0){
162// cout << "we are fitting " << string_1.first << " " << string_2.first << endl;
163 h1->Fit(&f1,"LQ", "same");
164
165// maxarr[ string.getIndex(string_1.first) ][ string.getIndex(string_2.first) ][0] = h1->GetXaxis()->GetBinCenter(h1->GetMaximumBin());
166// maxarr[ string.getIndex(string_1.first) ][ string.getIndex(string_2.first) ][1] = h1->GetMaximum();
167 maxarr[ string.getIndex(string_1.first) ][ string.getIndex(string_2.first) ][0] = h1->GetFunction("f1")->GetParameter(1);
168 maxarr[ string.getIndex(string_1.first) ][ string.getIndex(string_2.first) ][1] = h1->GetFunction("f1")->GetParameter(0);
169 maxarr[ string.getIndex(string_1.first) ][ string.getIndex(string_2.first) ][2] = h1->GetFunction("f1")->GetParameter(2);
170 }
171 }
172 } // if string 1 != string 2
173 } // loop over string 2
174 } // loop over modules of string 1
175 } // loop over string 1
176
177
178 for (int i=0;i<number_of_strings;i++){
179 for (int j=i+1;j<number_of_strings;j++){
180 if (maxarr[i][j][1]>setmax && maxarr[j][i][1]>setmax && maxarr[i][j][2]>20. && maxarr[j][i][2]>20.){
181 cout << i << " " << j << " " << maxarr[i][j][0] << " " << maxarr[j][i][0] << " " << (maxarr[i][j][0]-maxarr[j][i][0])/2. << " " << maxarr[i][j][1] << " " << maxarr[j][i][1] << " " << " " << duarr[i] << " " << duarr[j] << " " << neighbour << endl;
182 }
183 }
184 }
185
186
187
188 TFile out(outputFile.c_str(), "recreate");
189
190 out.cd();
191
192 for (auto& h1 : H1) {
193 for (auto& h2 : h1.second) {
194 h2.second.p->Write();
195 }
196 }
197
198 out.Close();
199}
string outputFile
Data structure for detector geometry and calibration.
#define FATAL(A)
Definition JMessage.hh:67
int debug
debug level
Definition JSirene.cc:72
int main(int argc, char **argv)
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition JParser.hh:2142
I/O formatting auxiliaries.
#define MAKE_CSTRING(A)
Make C-string.
Definition JPrint.hh:72
Direct access to string in detector data structure.
Detector data structure.
Definition JDetector.hh:96
General exception.
Definition JException.hh:24
Utility class to parse command line options.
Definition JParser.hh:1698
int read(const int argc, const char *const argv[])
Parse the program's command line options.
Definition JParser.hh:1992
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Detector file.
Definition JHead.hh:227
Router for mapping of string identifier to index.