Jpp 20.0.0-rc.3
the software that should make you happy
Loading...
Searching...
No Matches
Classes | Namespaces | Typedefs | Functions
JBillabong.cc File Reference

Main program to trigger acoustic data. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include <algorithm>
#include <chrono>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JModuleSupportkit.hh"
#include "JLang/JPredicate.hh"
#include "JLang/JComparator.hh"
#include "JTools/JHashMap.hh"
#include "JTools/JAbstractHistogram.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMeta.hh"
#include "JGeometry3D/JCylinder3D.hh"
#include "JGeometry3D/JOmega3D.hh"
#include "JGeometry3D/JRotator3D.hh"
#include "JAcoustics/JToA.hh"
#include "JAcoustics/JTransmission.hh"
#include "JAcoustics/JReceiver.hh"
#include "JAcoustics/JSoundVelocity.hh"
#include "JAcoustics/JAcousticsToolkit.hh"
#include "JAcoustics/JEvent.hh"
#include "JAcoustics/JEventOverlap.hh"
#include "JAcoustics/JSupport.hh"
#include "Jeep/JProperties.hh"
#include "Jeep/JPrint.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"
#include "JTrigger/JMatch.hh"
#include "JTrigger/JAlgorithm.hh"
#include "JTrigger/JBind2nd.hh"

Go to the source code of this file.

Classes

struct  JACOUSTICS::hit_type
 Transmission with position. More...
 
class  JACOUSTICS::JMatch3D
 3D match criterion for acoustic signals. More...
 
class  JACOUSTICS::JMatch1D
 1D match criterion for acoustic signals. More...
 
struct  JACOUSTICS::JCheck
 Auxiliary data structure for final check on event. More...
 

Namespaces

namespace  JACOUSTICS
 Auxiliary classes and methods for acoustic position calibration.
 

Typedefs

typedef JAbstractHistogram< double > JACOUSTICS::JHistogram_t
 Type definition for scan along axis.
 

Functions

void JACOUSTICS::print (std::ostream &out, const JEvent &event)
 Print event.
 
int main (int argc, char **argv)
 

Detailed Description

Main program to trigger acoustic data.

An acoustic event is based on coincidences between times of arrival.
If the number of coincident times of arrival exceeds a preset minimum, the event is triggered and subsequently stored in the output file.

Author
mdejong

Definition in file JBillabong.cc.

Function Documentation

◆ main()

int main ( int argc,
char ** argv )

Definition at line 315 of file JBillabong.cc.

316{
317 using namespace std;
318 using namespace JPP;
319
320 typedef JTYPELIST<JEvent, TH1D, JMeta>::typelist typelist;
321
322 JMultipleFileScanner<JToA> inputFile;
323 JLimit_t& numberOfEvents = inputFile.getLimit();
324
325 struct {
326 int factoryLimit = 100000;
327 int numberOfHits = 0;
328 double Q = 0.0;
329 double TMax_s = 0.0;
330 double DMax_m = 500.0;
331 double ZMax_m = 200.0;
332 double TMaxExtra_s = 100.0e-6;
333 double TMaxVertex_s = 1.0e-2;
334 double gridAngle_deg = 10.0;
335 } parameters;
336
337 JHistogram_t X;
338 JHistogram_t Y;
339 JHistogram_t Z;
340 JFileRecorder<typelist> outputFile;
341 JSoundVelocity V = getSoundVelocity; // default sound velocity
342 string detectorFile;
343 int waveform;
344 int debug;
345
346 try {
347
348 JProperties properties;
349
350 properties.insert(gmake_property(parameters.factoryLimit));
351 properties.insert(gmake_property(parameters.numberOfHits));
352 properties.insert(gmake_property(parameters.Q));
353 properties.insert(gmake_property(parameters.TMax_s));
354 properties.insert(gmake_property(parameters.DMax_m));
355 properties.insert(gmake_property(parameters.ZMax_m));
356 properties.insert(gmake_property(parameters.TMaxExtra_s));
357 properties.insert(gmake_property(parameters.TMaxVertex_s));
358 properties.insert(gmake_property(parameters.gridAngle_deg));
359
360 JParser<> zap("Main program to trigger acoustic data.");
361
362 zap['f'] = make_field(inputFile);
363 zap['n'] = make_field(numberOfEvents) = JLimit::max();
364 zap['@'] = make_field(properties, "trigger parameters") = JPARSER::initialised();
365 zap['X'] = make_field(X, "x-scan");
366 zap['Y'] = make_field(Y, "y-scan");
367 zap['Z'] = make_field(Z, "z-scan");
368 zap['o'] = make_field(outputFile, "output file") = "event.root";
369 zap['V'] = make_field(V, "sound velocity") = JPARSER::initialised();
370 zap['a'] = make_field(detectorFile, "detector file");
371 zap['W'] = make_field(waveform, "waveform identifier");
372 zap['d'] = make_field(debug) = 1;
373
374 zap(argc, argv);
375 }
376 catch(const exception &error) {
377 FATAL(error.what() << endl);
378 }
379
380 JDetector detector;
381
382 try {
383 load(detectorFile, detector);
384 }
385 catch(const JException& error) {
386 FATAL(error);
387 }
388
389 V.set(detector.getUTMZ());
390
391 const JCylinder3D cylinder(detector.begin(), detector.end());
392
393 const double v0 = V(cylinder.getZmax());
394
395 if (parameters.TMax_s < parameters.DMax_m / v0) {
396
397 parameters.TMax_s = parameters.DMax_m / v0;
398
399 NOTICE("Set maximal time [s] " << FIXED(7,3) << parameters.TMax_s << endl);
400 }
401
402 const JMatch3D match3D(v0, parameters.DMax_m, parameters.TMaxExtra_s);
403 const JMatch1D match1D(v0, parameters.DMax_m, parameters.ZMax_m, parameters.TMaxExtra_s);
404
405 const JEventOverlap overlap(parameters.TMax_s);
406
407 const JCheck check(X, Y, Z, v0, parameters.TMaxVertex_s, parameters.numberOfHits);
408
409 const JOmega3D omega(JAngle3D(0.0, 0.0), JOmega3D::range_type(0.0, 0.5*PI), parameters.gridAngle_deg * PI / 180.0);
410 const JRotator3D rotator(omega);
411
412
413 JHashMap<int, JReceiver> receivers;
414
415 const JPosition3D center(cylinder.getX(), cylinder.getY(), 0.5 * (cylinder.getZmin() + cylinder.getZmax()));
416
417 for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
418
419 if (module->getFloor() != 0) {
420
421 receivers[module->getID()] = JReceiver(module->getID(),
422 module->getPosition() + getPiezoPosition() - center,
423 module->getT0() * 1.0e-9);
424 }
425 }
426
427 TH1D h0("h0", NULL, 21, -0.5, 20.5);
428
429 vector<TH1D*> M1(10, NULL);
430
431 for (size_t i = 2; i <= 5; ++i) {
432 M1[i] = new TH1D(MAKE_CSTRING("M[" << i << "]"), NULL, 1001, -0.5, 1000.5);
433 }
434
435 outputFile.open();
436
437 if (!outputFile.is_open()) {
438 FATAL("Error opening file " << outputFile << endl);
439 }
440
441 outputFile.put(JMeta(argc, argv));
442
443 // input data
444
445 typedef vector<hit_type> buffer_type; // data type
446
447 buffer_type data;
448
449 while (inputFile.hasNext()) {
450
451 if (inputFile.getCounter()%1000 == 0) {
452 STATUS("counter: " << setw(8) << inputFile.getCounter() << '\r' << flush); DEBUG(endl);
453 }
454
455 const JToA* p = inputFile.next();
456
457 if (detector.getID() != p->DETID) { // consistency check
458 FATAL("Invalid detector identifier " << p->DETID << " != " << detector.getID() << endl);
459 }
460
461 if (p->WAVEFORMID == waveform && receivers.has(p->DOMID)) {
462
463 if (p->QUALITYFACTOR >= parameters.Q * (parameters.Q <= 1.0 ? p->QUALITYNORMALISATION : 1.0)) {
464
465 const JReceiver& receiver = receivers[p->DOMID];
466
467 const double toa = p->TOA_S();
468
469 const JTransmission transmission(p->RUN,
470 p->DOMID,
471 p->QUALITYFACTOR,
472 p->QUALITYNORMALISATION,
473 receiver.getT(toa),
474 receiver.getT(toa));
475
476 data.push_back(hit_type(receiver.getPosition(), transmission));
477 }
478 }
479 }
480 STATUS(endl);
481
482
483 sort(data.begin(), data.end(), make_comparator(&hit_type::getToA, JComparison::lt())); // sort according time-of-arrival
484
485
486 buffer_type buffer(parameters.factoryLimit); // local buffer of hits
487 JEvent out[2]; // FIFO of events
488
489 const chrono::steady_clock::time_point t0 = chrono::steady_clock::now();
490
491 size_t trigger_counter = 0;
492
493 for (buffer_type::const_iterator p = data.begin(), q = p; p != data.end(); ++p) {
494
495 if (distance(data.cbegin(),p)%10000 == 0) {
496 //STATUS("processed: " << FIXED(5,1) << (double) (100 * distance(data.cbegin(),p)) / (double) data.size() << "%" << ' ' << FIXED(12,2) << p->getToA() << " [s]" << '\r' << flush); DEBUG(endl);
497
498 const chrono::steady_clock::time_point t1 = chrono::steady_clock::now();
499
500 STATUS("processed: "
501 << FIXED(5,1) << (double) (100 * distance(data.cbegin(),p)) / (double) data.size() << "%" << ' '
502 << FIXED(8,1) << p->getToA() << " [s]" << ' '
503 << setw(6) << chrono::duration_cast<chrono::seconds>(t1 - t0).count() << " [s]" << ' '
504 << setw(6) << trigger_counter << endl);
505 }
506
507 h0.Fill(0.0);
508
509 while (q != data.end() && q->getToA() - p->getToA() <= parameters.TMax_s) { ++q; }
510
511 if (distance(p,q) >= parameters.numberOfHits) {
512
513 h0.Fill(1.0);
514
515 if (distance(p,q) < parameters.factoryLimit) {
516
517 h0.Fill(2.0);
518
519 buffer_type::iterator root = buffer.begin();
520 buffer_type::iterator __p = buffer.begin();
521 buffer_type::iterator __q = buffer.begin();
522
523 *root = *p;
524
525 ++__p;
526 ++__q;
527
528 for (buffer_type::const_iterator i = p; ++i != q; ) {
529 if (match3D(*root,*i)) {
530 *(__q++) = *i;
531 }
532 }
533
534 M1[2]->Fill(distance(root,__q));
535
536 if (distance(root,__q) >= parameters.numberOfHits) {
537
538 h0.Fill(3.0);
539
540 __q = clusterize(__p, __q, match3D);
541
542 M1[3]->Fill(distance(root,__q));
543
544 if (distance(root,__q) >= parameters.numberOfHits) {
545
546 h0.Fill(4.0);
547
548 const double W = 1.0 / (double) rotator.size();
549
550 for (const JRotation3D& R : rotator) {
551
552 for (buffer_type::iterator i = root; i != __q; ++i) {
553 i->rotate(R);
554 }
555
556 buffer_type::iterator __z = partition(__p, __q, JBind2nd(match1D,*root));
557
558 M1[4]->Fill(distance(root,__z), W);
559
560 if (distance(root,__z) >= parameters.numberOfHits) {
561
562 h0.Fill(5.0, W);
563
564 __z = clusterize(__p, __z, match1D);
565
566 M1[5]->Fill(distance(root,__z), W);
567
568 if (distance(root,__z) >= parameters.numberOfHits) {
569
570 h0.Fill(6.0, W);
571
572 if (!check.is_valid() || check(*root, __p, __z)) {
573
574 trigger_counter += 1;
575
576 h0.Fill(7.0, W);
577
578 out[1] = JEvent(detector.getID(), out[0].getCounter() + 1, waveform, p, q);
579
580 break;
581 }
582 }
583 }
584 }
585 }
586 }
587
588 } else {
589
590 trigger_counter += 1;
591
592 out[1] = JEvent(detector.getID(), out[0].getCounter() + 1, waveform, p, q);
593 }
594
595 if (!out[1].empty()) {
596
597 if (out[0].empty()) {
598
599 out[0] = out[1]; // shift
600
601 } else if (overlap(out[0],out[1])) {
602
603 out[0].merge(out[1]); // merge
604
605 } else {
606
607 if (debug >= notice_t) { print(cout, out[0]); }
608
609 outputFile.put(out[0]); // write
610
611 out[0] = out[1]; // shift
612 }
613 }
614
615 out[1].clear();
616 }
617 }
618
619 if (!out[0].empty()) {
620
621 if (debug >= notice_t) { print(cout, out[0]); }
622
623 outputFile.put(out[0]); // write
624 }
625 STATUS(endl);
626
627 outputFile.put(h0);
628
629 for (TH1D* h1 : M1) {
630 if (h1 != NULL) {
631 outputFile.put(*h1);
632 }
633 }
634
635 JMultipleFileScanner<JMeta> io(inputFile);
636
637 io >> outputFile;
638
639 outputFile.close();
640}
outputFile
string outputFile
Definition JDAQTimesliceSelector.cc:37
DEBUG
#define DEBUG(A)
Message macros.
Definition JMessage.hh:62
STATUS
#define STATUS(A)
Definition JMessage.hh:63
NOTICE
#define NOTICE(A)
Definition JMessage.hh:64
FATAL
#define FATAL(A)
Definition JMessage.hh:67
debug
int debug
debug level
Definition JSirene.cc:72
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition JParser.hh:2142
MAKE_CSTRING
#define MAKE_CSTRING(A)
Make C-string.
Definition JPrint.hh:72
gmake_property
#define gmake_property(A)
macros to convert (template) parameter to JPropertiesElement object
Definition JProperties.hh:1139
distance
std::vector< T >::difference_type distance(typename std::vector< T >::const_iterator first, typename PhysicsEvent::const_iterator< T > second)
Specialisation of STL distance.
Definition build/antares-dataformat/PhysicsEvent.hh:434
JACOUSTICS::JMatch3D
3D match criterion for acoustic signals.
Definition JBillabong.cc:96
JDETECTOR::JDetector
Detector data structure.
Definition JDetector.hh:96
JEEP::JProperties
Utility class to parse parameter values.
Definition JProperties.hh:503
JGEOMETRY3D::JAngle3D
Data structure for angles in three dimensions.
Definition JAngle3D.hh:35
JGEOMETRY3D::JOmega3D
Direction set covering (part of) solid angle.
Definition JOmega3D.hh:68
JGEOMETRY3D::JPosition3D
Data structure for position in three dimensions.
Definition JPosition3D.hh:38
JGEOMETRY3D::JPosition3D::getPosition
const JPosition3D & getPosition() const
Get position.
Definition JPosition3D.hh:130
JLANG::JException
General exception.
Definition JException.hh:24
JPARSER::JParser
Utility class to parse command line options.
Definition JParser.hh:1698
JSUPPORT::JFileRecorder
Object writing to file.
Definition JFileRecorder.hh:44
JSUPPORT::JMultipleFileScanner
General purpose class for object reading from a list of file names.
Definition JMultipleFileScanner.hh:306
JSUPPORT::JMultipleFileScanner::hasNext
virtual bool hasNext() override
Check availability of next element.
Definition JMultipleFileScanner.hh:409
JSUPPORT::JMultipleFileScanner::getCounter
counter_type getCounter() const
Get counter.
Definition JMultipleFileScanner.hh:382
JSUPPORT::JMultipleFileScanner::next
virtual const pointer_type & next() override
Get next element.
Definition JMultipleFileScanner.hh:457
JTOOLS::JHashCollection::has
bool has(const T &value) const
Test whether given value is present.
Definition JHashCollection.hh:309
JTRIGGER::JMatch1D
1D match criterion.
Definition JMatch1D.hh:33
JACOUSTICS::getSoundVelocity
static const JSoundVelocity getSoundVelocity(1541.0, -17.0e-3, -2000.0)
Function object for velocity of sound.
JCOMPAREHISTOGRAMS::print
std::ostream & print(std::ostream &out, const JTestSummary &summary, T __begin, T __end, const bool useColors=true, const JFormat_t &formatting=JFormat_t(18, 3, std::ios::fixed))
Print test summary.
Definition JTestSummary.hh:254
JDETECTOR::load
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Definition JDetectorToolkit.hh:536
JDETECTOR::getPiezoPosition
JPosition3D getPiezoPosition()
Get relative position of piezo in optical module.
Definition JModuleSupportkit.hh:30
JEEP::notice_t
@ notice_t
notice
Definition JMessage.hh:32
JLANG::make_comparator
JComparator< JResult_t T::*, JComparison::lt > make_comparator(JResult_t T::*member)
Helper method to create comparator between values of data member.
Definition JLang/JComparator.hh:185
JPP
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Definition JAAnetToolkit.hh:45
JRECONSTRUCTION::buffer_type
std::vector< JHitW0 > buffer_type
hits
Definition JPerth.cc:70
JRECONSTRUCTION::JEvent
JFIT::JEvent JEvent
Definition JHistory.hh:454
JTOOLS::overlap
bool overlap(const JRange< T, JComparator_t > &first, const JRange< T, JComparator_t > &second)
Test overlap between ranges.
Definition JRange.hh:641
JTRIGGER::JBind2nd
JBinder2nd< JHit_t > JBind2nd(const JMatch< JHit_t > &match, const JHit_t &second)
Auxiliary method to create JBinder2nd object.
Definition JBind2nd.hh:66
JTRIGGER::clusterize
JHitIterator_t clusterize(JHitIterator_t __begin, JHitIterator_t __end, const JMatch_t &match, const int Nmin=1)
Partition data according given binary match operator.
Definition JAlgorithm.hh:38
root
Definition root.py:1
FIXED
Auxiliary data structure for floating point format specification.
Definition JManip.hh:448
JAANET::detector
Detector file.
Definition JHead.hh:227
JACOUSTICS::JCheck
Auxiliary data structure for final check on event.
Definition JBillabong.cc:200
JACOUSTICS::JEventOverlap
Match of two events considering overlap in time and position.
Definition JAcoustics/JEventOverlap.hh:20
JACOUSTICS::JEvent::merge
void merge(const JEvent &event)
Merge event.
Definition JAcoustics/JEvent.hh:165
JACOUSTICS::JReceiver
Acoustic receiver.
Definition JReceiver.hh:30
JACOUSTICS::JReceiver::getT
double getT(const double t_s) const
Get corrected time.
Definition JReceiver.hh:72
JACOUSTICS::JSoundVelocity
Implementation for depth dependend velocity of sound.
Definition JSoundVelocity.hh:37
JACOUSTICS::JSoundVelocity::set
JSoundVelocity & set(const double z0)
Set depth.
Definition JSoundVelocity.hh:62
JACOUSTICS::JToA
Time-of-arrival data from acoustic piezo sensor or hydrophone.
Definition JToA.hh:28
JACOUSTICS::JToA::DOMID
uint32_t DOMID
DAQ run number.
Definition JToA.hh:34
JACOUSTICS::JToA::QUALITYFACTOR
uint32_t QUALITYFACTOR
The ticks (16ns) part of the DAQ frame timestamp.
Definition JToA.hh:39
JACOUSTICS::JToA::QUALITYNORMALISATION
uint32_t QUALITYNORMALISATION
A measure of how good the waveform match was to the signal.
Definition JToA.hh:40
JACOUSTICS::JToA::WAVEFORMID
int32_t WAVEFORMID
DOM unique identifeir.
Definition JToA.hh:35
JACOUSTICS::JToA::DETID
int32_t DETID
Definition JToA.hh:32
JACOUSTICS::JToA::RUN
int32_t RUN
detector identifier
Definition JToA.hh:33
JACOUSTICS::JToA::TOA_S
double TOA_S() const
Time of Arrival, expressed in seconds relative to Unix epoch (1 January 1970 00:00:00 UTC)
Definition JToAImp.cc:40
JACOUSTICS::JTransmission
Acoustic transmission.
Definition JTransmission.hh:31
JACOUSTICS::hit_type
Transmission with position.
Definition JBillabong.cc:70
JLANG::JTypeList
Type list.
Definition JTypeList.hh:23
JPARSER::initialised
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition JParser.hh:68
JSUPPORT::JLimit
Auxiliary class for defining the range of iterations of objects.
Definition JLimit.hh:45
JSUPPORT::JLimit::max
static counter_type max()
Get maximum counter value.
Definition JLimit.hh:128
JSUPPORT::JMeta
Auxiliary class for ROOT I/O of application specific meta data.
Definition JMeta.hh:72
JTOOLS::JAbstractHistogram
Simple data structure for histogram binning.
Definition JAbstractHistogram.hh:26
JTOOLS::JHashMap
General purpose class for hash map of unique keys.
Definition JHashMap.hh:75
Generated by doxygen 1.12.0