Jpp 19.3.0-rc.3
the software that should make you happy
Loading...
Searching...
No Matches
JDepthDependence.cc
Go to the documentation of this file.
1#include <iomanip>
2#include <iostream>
3#include <sstream>
4#include <string>
5
6#include "TFile.h"
7#include "TH1D.h"
8#include "TMath.h"
9#include "TROOT.h"
10
11#include "TGraphErrors.h"
12
17#include "Jeep/JParser.hh"
18#include "Jeep/JMessage.hh"
19#include "JLang/JString.hh"
21#include "JSupport/JSupport.hh"
22
23
24namespace {
25
26 using namespace JPP;
27 using namespace std;
28
29 double getCount (TH1D* hptr, int muon_threshold) {
30 double count = 0;
31 for (int i = hptr->GetXaxis()->FindBin(muon_threshold); i <= hptr->GetNbinsX(); i++) {
32 count += hptr->GetBinContent(i);
33 }
34 return count;
35 }
36
37 double getLiveTime (TH1D* hptr, const JModule& module) {
38 string label = getLabel(module);
39 double livetime = hptr->GetBinContent(hptr->GetXaxis()->FindBin(label.c_str()));
40 return livetime;
41 }
42
43}
44
45/**
46 * \file
47 *
48 * Auxiliary program to find depth dependence of multiplicity rates
49 * \author mlincetto
50 */
51int main(int argc, char **argv)
52{
53 using namespace std;
54 using namespace JPP;
55
56 //-----------------------------------------------------------
57 // parameter interface
58 //-----------------------------------------------------------
59
60
61 string inputFile;
62 string outputFile;
63 string detectorFile;
64 int line;
65 int debug;
66 int minMultiplicity;
67 int minLiveTime_s;
68 // string summaryFile;
69
70 try {
71
72 JParser<> zap("Auxiliary program to find depth dependence of multiplicity rates");
73
74 zap['f'] = make_field(inputFile, "JMM input file");
75 zap['o'] = make_field(outputFile, "output file") = "depthdependence.root";
76 // zap['s'] = make_field(summaryFile, "summary file") = "depthdependence.txt";
77 zap['a'] = make_field(detectorFile);
78 zap['d'] = make_field(debug) = 1;
79 zap['l'] = make_field(line) = 1;
80 zap['M'] = make_field(minMultiplicity, "Minimum multiplicity") = 8;
81 zap['T'] = make_field(minLiveTime_s, "Minimum DOM livetime [s] to be eligible for plotting") = 3600;
82
83 zap(argc, argv);
84 }
85
86 catch(const exception &error) {
87 FATAL(error.what() << endl);
88 }
89
90
91 //----------------------
92 // loading detector file
93 //----------------------
94
96
97 try {
98 load(detectorFile, detector);
99 }
100 catch(const JException& error) {
101 FATAL(error);
102 }
103
104 if (detector.empty())
105 FATAL("Empty detector." << endl);
106
107 const JModuleRouter router(detector);
108
109 double utm_z = detector.getUTMZ();
110
111 NOTICE("Detector base UTM z [m]: " << utm_z << endl);
112
113 //-----------------------------------------------------------
114 // load input files
115 //-----------------------------------------------------------
116
117
118 DEBUG("Loading input file " << inputFile << endl);
119
120 TFile* in = TFile::Open(inputFile.c_str(), "exist");
121
122 if (in == NULL || !in->IsOpen()) {
123 FATAL("File: " << inputFile << " not opened." << endl);
124 }
125
126 //-----------------------------------------------------------
127 // recover rates from histograms
128 //-----------------------------------------------------------
129
130
131 DEBUG("Loading livetime histogram from " << inputFile << endl);
132
133 TH1D* liveTime = (TH1D*)in->Get("LiveTime");
134
135 if (liveTime == NULL) {
136 FATAL("Missing live time histogram.");
137 }
138
140 vector<double> rate_val;
141 vector<double> rate_err;
142
143 for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
144
145 if (module->getString() != line) { continue; }
146
147 STATUS("Loading data histogram for from input file: " << getLabel(*module) << "\t ID = " << module->getID() << endl);
148
149 TH1D* data_histogram = (TH1D*)in->Get(TString(getLabel(*module)) + "_P");
150
151 if (data_histogram != NULL) {
152
153 double data_count = getCount(data_histogram, minMultiplicity);
154 double data_livetime = getLiveTime(liveTime, *module);
155
156 double module_depth = - utm_z - module->getZ();
157
158 if (data_livetime > minLiveTime_s) {
159
160 double val = data_count / data_livetime;
161 double err = sqrt(data_count) / data_livetime;
162
163 NOTICE(module_depth << " " << val << endl);
164
165 rate_val.push_back(val);
166 rate_err.push_back(err);
167 depth.push_back( module_depth);
168
169 }
170
171 }
172
173 }
174
175
176 TGraph* gr_data = new TGraphErrors(depth.size(), &depth[0], &rate_val[0], 0, &rate_err[0]);
177
178 gr_data->SetTitle(TString("KM3NeT Preliminary; Depth [m]; Inclusive ") + Form("%d", minMultiplicity) + TString("-fold coincidence rate [Hz]"));
179
180 TFile out(outputFile.c_str(), "recreate");
181
182 gr_data->Write();
183
184 out.Close();
185
186 return 0;
187
188}
string outputFile
KM3NeT DAQ constants, bit handling, etc.
int main(int argc, char **argv)
Data structure for detector geometry and calibration.
General purpose messaging.
#define DEBUG(A)
Message macros.
Definition JMessage.hh:62
#define STATUS(A)
Definition JMessage.hh:63
#define NOTICE(A)
Definition JMessage.hh:64
#define FATAL(A)
Definition JMessage.hh:67
int debug
debug level
Definition JSirene.cc:72
Direct access to module in detector data structure.
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Utility class to parse command line options.
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition JParser.hh:2142
ROOT TTree parameter settings of various packages.
Detector data structure.
Definition JDetector.hh:96
Router for direct addressing of module data in detector data structure.
Data structure for a composite optical module.
Definition JModule.hh:75
General exception.
Definition JException.hh:24
Utility class to parse command line options.
Definition JParser.hh:1698
std::string getLabel(const JLocation &location)
Get module label for monitoring and other applications.
Definition JLocation.hh:247
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
size_t getCount(const array_type< T > &buffer, const JCompare_t &compare)
Count number of unique values.
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Depth.
Definition JHead.hh:980
Detector file.
Definition JHead.hh:227
Normalisation of MUPAGE events.
Definition JHead.hh:835