13 #include "evt/Head.hh"
72 inline bool compare(
const JHitL0& first,
const JHitL0& second)
74 if (std::equal_to<KM3NETDAQ::JDAQPMTIdentifier>()(first, second))
75 return std::less<JTRIGGER::JHit>()(first, second);
77 return std::less<KM3NETDAQ::JDAQPMTIdentifier>()(first, second);
88 int main(
int argc,
char **argv)
92 using namespace KM3NETDAQ;
95 typedef JParallelFileScanner_t::multi_pointer_type multi_pointer_type;
98 JParallelFileScanner_t inputFile;
105 size_t numberOfPrefits;
114 JParser<> zap(
"Program to perform PDF fit of muon trajectory to data.");
119 zap[
'n'] =
make_field(numberOfEvents) = JLimit::max();
121 zap[
'R'] =
make_field(roadWidth_m) = numeric_limits<double>::max();
127 zap[
'c'] =
make_field(calibrate,
"Enable time calibration per optical module.");
132 catch(
const exception& error) {
133 FATAL(error.what() << endl);
161 JRegressor_t::T_ns.setRange(-50.0, +450.0);
162 JRegressor_t::Vmax_npe = 10.0;
163 JRegressor_t::MAXIMUM_ITERATIONS = 1000;
165 JRegressor_t fit(pdfFile, TTS_ns);
167 fit.parameters.resize(5);
169 fit.parameters[0] = JLine3Z::pX();
170 fit.parameters[1] = JLine3Z::pY();
171 fit.parameters[2] = JLine3Z::pT();
172 fit.parameters[3] = JLine3Z::pDX();
173 fit.parameters[4] = JLine3Z::pDY();
175 if (fit.getRmax() < roadWidth_m) {
176 roadWidth_m = fit.getRmax();
179 TH2D h2(
"h2", NULL,
detector.size(), -0.5,
detector.size() - 0.5, 1001, -250.0, +250.0);
187 while (inputFile.hasNext()) {
189 STATUS(
"event: " << setw(10) << inputFile.getCounter() <<
'\r');
DEBUG(endl);
191 multi_pointer_type ps = inputFile.next();
197 JEvt::iterator __end = evt->end();
203 if (evt->begin() != __end) {
205 copy(evt->begin(), __end, back_inserter(out));
207 if (numberOfPrefits > 0) {
208 advance(__end = evt->begin(), min(numberOfPrefits, out.size()));
211 partial_sort(evt->begin(), __end, evt->end(),
qualitySorter);
213 __end = partition(evt->begin(), __end,
JHistory::is_event(evt->begin()->getHistory()));
218 buildL0(*tev, moduleRouter,
true, back_inserter(dataL0));
221 if (dataL0.size() >= fit.parameters.size()) {
223 for (JEvt::const_iterator track = evt->begin(); track != __end; ++track) {
233 for (JDataL0_t::const_iterator i = dataL0.begin(); i != dataL0.end(); ++i) {
246 sort(data.begin(), data.end(), compare);
248 JDataW0_t::iterator __end = unique(data.begin(), data.end(), equal_to<JDAQPMTIdentifier>());
251 const int NDF = distance(data.begin(), __end) - fit.parameters.size();
257 if (track->getE() > 0.1)
258 fit.E_GeV = track->getE();
262 const double chi2 = fit(
JLine3Z(tz), data.begin(), __end);
273 fit.error.getDY() * fit.error.getDY()));
285 for (JDataW0_t::const_iterator hit = data.begin(); hit != __end; ++hit) {
286 buffer.insert(hit->getModuleID());
295 if (distance(data.begin(), q) - fit.parameters.size() > 0) {
297 fit(ta, data.begin(), q);
299 sort(q, __end, compare);
301 h2.Fill(moduleRouter.
getIndex(*
id), q->getT() - fit.value.getT(q->getPosition()));
Data regression method for JFIT::JLine3Z.
Utility class to parse command line options.
General purpose data regression method.
JPredicate< JResult_t T::*, JComparison::eq > make_predicate(JResult_t T::*member, const JResult_t value)
Helper method to create predicate for data member.
number of iterations from JGandalf.cc
Template specialisation of L0 builder for JHitL0 data type.
const int getIndex(const JObjectID &id) const
Get index of module.
Recording of objects on file according a format that follows from the file name extension.
Router for direct addressing of module data in detector data structure.
General purpose class for parallel reading of objects from a single file or multiple files...
Template specialisation of class JModel to match hit with muon trajectory along z-axis.
Auxiliary class to test history.
Data structure for fit of straight line in positive z-direction.
Container for historical events.
Data structure for detector geometry and calibration.
angular resolution [rad] from JGandalf.cc
Various implementations of functional maps.
JFit getFit(const JHistory &history, const JTrack3D &track, const double Q, const int NDF, const double energy=0.0, const int status=0)
Get fit.
Basic data structure for L0 hit.
Definition of fit parameters from various applications.
Basic data structure for time and time over threshold information of hit.
Auxiliary class for defining the range of iterations of objects.
JAxis3D & rotate(const JRotation3D &R)
Rotate axis.
control parameter from JGandalf.cc
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Auxiliary class to test history.
Parallel scanning of objects from a single file or multiple files according a format that follows fro...
Auxiliary class for a hit with background rate value.
angular resolution [rad] from JGandalf.cc
void load(const JString &file_name, JDetector &detector)
Load detector from input file.
General purpose messaging.
counter_type advance(counter_type &counter, const counter_type value, const counter_type limit=std::numeric_limits< counter_type >::max())
Advance counter.
Scanning of objects from multiple files according a format that follows from the extension of each fi...
Direct access to module in detector data structure.
number of hits from JGandalf.cc
double getQuality(const double chi2, const int NDF)
Get quality of fit.
Data structure for set of track fit results.
Regressor function object for JLine3Z fit using JGandalf minimiser.
Utility class to parse command line options.
Data structure for L0 hit.
ROOT TTree parameter settings.
Data structure for fit of straight line paralel to z-axis.
bool qualitySorter(const JFIT::JFit &first, const JFIT::JFit &second)
Comparison of fit results.
void copy(const Head &from, JHead &to)
Copy header from from to to.
JDirection3D getDirection(const Vec &v)
Get direction.
Object reading from a list of files.
const JLimit & getLimit() const
Get limit.
JPosition3D & rotate(const JRotation3D &R)
Rotate.
#define DEBUG(A)
Message macros.
JPosition3D getPosition(const Vec &v)
Get position.
int main(int argc, char *argv[])