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

Main program to trigger acoustic data. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include <set>
#include <map>
#include <algorithm>
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "TGraph.h"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JHydrophone.hh"
#include "JLang/JClonable.hh"
#include "JLang/JPredicate.hh"
#include "JLang/JComparator.hh"
#include "JTools/JHashMap.hh"
#include "JROOT/JManager.hh"
#include "JROOT/JGraph.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMeta.hh"
#include "JAcoustics/JToA.hh"
#include "JAcoustics/JTransmission.hh"
#include "JAcoustics/JReceiver.hh"
#include "JAcoustics/JSoundVelocity.hh"
#include "JAcoustics/JConstantSoundVelocity.hh"
#include "JAcoustics/JAcousticsToolkit.hh"
#include "JAcoustics/JAcousticsSupportkit.hh"
#include "JAcoustics/JTriggerParameters.hh"
#include "JAcoustics/JEvent.hh"
#include "JAcoustics/JSupport.hh"
#include "JGeometry3D/JVector3D.hh"
#include "JGeometry3D/JVertex3D.hh"
#include "JGeometry3D/JPosition3D.hh"
#include "JGeometry3D/JCylinder3D.hh"
#include "Jeep/JContainer.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"

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...
 
struct  JACOUSTICS::event_type
 Event with vertex. More...
 
struct  JACOUSTICS::JEventOverlap
 Match of two events considering overlap in time and position. More...
 
struct  JACOUSTICS::vertex_type
 Vertex. More...
 
struct  JACOUSTICS::JVelo
 Vertex locator. More...
 

Namespaces

namespace  JACOUSTICS
 Auxiliary classes and methods for acoustic position calibration.
 

Functions

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 JAcousticsTriggerProcessor.cc.

Function Documentation

◆ main()

int main ( int argc,
char ** argv )

Definition at line 423 of file JAcousticsTriggerProcessor.cc.

424{
425 using namespace std;
426 using namespace JPP;
427
428 typedef JTYPELIST<JEvent, JTriggerParameters, JMeta, TH1D, TGraph>::typelist typelist;
429
430 typedef JContainer< vector<JHydrophone> > hydrophones_container;
431
432 JMultipleFileScanner<JToA> inputFile;
433 JLimit_t& numberOfEvents = inputFile.getLimit();
434 JTriggerParameters parameters;
435 int factoryLimit = 10000;
436 double TMaxExtra_s = 500.0e-6;
437 double RMax_m = 0.0;
438 double DMax_m = 0.0;
439 double Xv_m = 0.0;
440 double factor = 2.0;
441 double Z_m = 0.0; // z-position of emitters
442 double fraction = 0.75; // fraction of working modules
443 double grid = 10.0;
444 JFileRecorder<typelist> outputFile;
445 JSoundVelocity V = getSoundVelocity; // default sound velocity
446 string detectorFile;
447 hydrophones_container hydrophones; // hydrophones
448 double precision;
449 set<int> ID;
450 int debug;
451
452 try {
453
454 JProperties properties;
455
456 properties.insert(gmake_property(parameters.Q));
457 properties.insert(gmake_property(parameters.TMax_s));
458 properties.insert(gmake_property(parameters.numberOfHits));
459 properties.insert(gmake_property(fraction));
460 properties.insert(gmake_property(factoryLimit));
461 properties.insert(gmake_property(TMaxExtra_s));
462 properties.insert(gmake_property(RMax_m));
463 properties.insert(gmake_property(DMax_m));
464 properties.insert(gmake_property(Xv_m));
465 properties.insert(gmake_property(factor));
466 properties.insert(gmake_property(Z_m));
467 properties.insert(gmake_property(grid));
468
469 JParser<> zap("Main program to trigger acoustic data.");
470
471 zap['f'] = make_field(inputFile, "output of JToA");
472 zap['n'] = make_field(numberOfEvents) = JLimit::max();
473 zap['@'] = make_field(properties, "trigger parameters");
474 zap['o'] = make_field(outputFile, "output file") = "event.root";
475 zap['V'] = make_field(V, "sound velocity") = JPARSER::initialised();
476 zap['a'] = make_field(detectorFile, "detector file");
477 zap['H'] = make_field(hydrophones, "hydrophone data") = JPARSER::initialised();
478 zap['p'] = make_field(precision, "precision time-of-arrival") = 1.0e-6;
479 zap['W'] = make_field(ID, "waveform identifiers") = JPARSER::initialised();
480 zap['d'] = make_field(debug) = 1;
481
482 zap(argc, argv);
483 }
484 catch(const exception &error) {
485 FATAL(error.what() << endl);
486 }
487
488
489 JDetector detector;
490
491 try {
492 load(detectorFile, detector);
493 }
494 catch(const JException& error) {
495 FATAL(error);
496 }
497
498 V.set(detector.getUTMZ());
499
500 const JCylinder3D cylinder(detector.begin(), detector.end());
501
502 const JVector3D center(cylinder.getX(), cylinder.getY(), Z_m);
503
504 //const JConstantSoundVelocity V0 = V(cylinder.getZmax());
505 const JAbstractSoundVelocity& V0 = V;
506
507 const JVelo Velo(V0, RMax_m, Xv_m*grid, factor);
508 const JVelo velo(V0, RMax_m/grid, Xv_m, factor);
509
510 const JMatch3D match(V0, TMaxExtra_s);
511 const JEventOverlap overlap2D(parameters.TMax_s, DMax_m);
512 const JEventOverlap overlap1D(parameters.TMax_s);
513
514 JHashMap<int, JReceiver> receivers;
515
516 for (JDetector::const_iterator i = detector.begin(); i != detector.end(); ++i) {
517
518 JPosition3D pos(0.0, 0.0, 0.0);
519
520 if (i->getFloor() == 0) { // get relative position of hydrophone
521
522 try {
523 pos = getPosition(hydrophones.begin(),
524 hydrophones.end(),
525 make_predicate(&JHydrophone::getString, i->getString()));
526 }
527 catch(const exception&) {
528 continue; // if no position available, discard hydrophone
529 }
530 }
531
532 receivers[i->getID()] = JReceiver(i->getID(),
533 i->getPosition() + pos,
534 i->getT0() * 1.0e-9);
535 }
536
537 // monitoring
538
539 JManager<int, TH1D> M1(new TH1D("M1[%]", NULL, 2001, -0.5, 2000.5));
540 JManager<int, TH1D> M2(new TH1D("M2[%]", NULL, 2001, -0.5, 2000.5));
541 JManager<int, TH1D> M3(new TH1D("M3[%]", NULL, 2001, -0.5, 2000.5));
542 JManager<int, TH1D> M4(new TH1D("M4[%]", NULL, 2001, -0.5, 2000.5));
543 JManager<int, TH1D> M5(new TH1D("M5[%]", NULL, 2001, -0.5, 2000.5));
544
545 map<int, JGraph_t> G2;
546
547 outputFile.open();
548
549 if (!outputFile.is_open()) {
550 FATAL("Error opening file " << outputFile << endl);
551 }
552
553 outputFile.put(JMeta(argc, argv));
554 outputFile.put(parameters);
555
556
557 // input data
558
559 typedef vector<hit_type> buffer_type; // data type
560
561 map<int, map<int, buffer_type> > f1; // emitter -> receiver -> data
562
563 set<int> ids; // working receivers
564
565 STATUS("reading data... " << flush); DEBUG(endl);
566
567 while (inputFile.hasNext()) {
568
569 JToA* parameters = inputFile.next();
570
571 if (detector.getID() != parameters->DETID) { // consistency check
572 FATAL("Invalid detector identifier " << parameters->DETID << " != " << detector.getID() << endl);
573 }
574
575 if (receivers.has(parameters->DOMID)) {
576
577 ids.insert(parameters->DOMID);
578
579 const JReceiver& receiver = receivers[parameters->DOMID];
580 const double toa_s = parameters->TOA_S();
581
582 if (toa_s < event_type::TOA_s) {
583 event_type::TOA_s = toa_s;
584 }
585
586 const JTransmission transmission(parameters->RUN,
587 parameters->DOMID,
588 parameters->QUALITYFACTOR,
589 parameters->QUALITYNORMALISATION,
590 receiver.getT(toa_s),
591 receiver.getT(toa_s));
592
593 if (ID.empty() || ID.count(getWaveformID(parameters->WAVEFORMID)) != 0) {
594 f1[getWaveformID(parameters->WAVEFORMID)][receiver.getID()].push_back(hit_type(receiver.getPosition(), transmission));
595 }
596 }
597 }
598 STATUS("OK" << endl);
599
600 if (parameters.numberOfHits == 0) {
601
602 parameters.numberOfHits = (int) (ids.size() * fraction);
603
604 STATUS("Number of hits " << parameters.numberOfHits << endl);
605 }
606
607 vector<event_type> queue;
608
609 for (map<int, map<int, buffer_type> >::iterator i = f1.begin(); i != f1.end(); ++i) {
610
611 STATUS("processing: " << setw(3) << i->first << endl);
612
613 buffer_type data;
614
615 for (map<int, buffer_type>::iterator receiver = i->second.begin(); receiver != i->second.end(); ++receiver) {
616
617 // filter similar hits
618
619 sort(receiver->second.begin(), receiver->second.end(), JTransmission::compare(precision));
620
621 buffer_type::iterator __end = unique(receiver->second.begin(), receiver->second.end(), JTransmission::equals(precision));
622
623 // selection based on quality
624
625 for (buffer_type::const_iterator p = receiver->second.begin(); p != __end; ++p) {
626 if (p->getQ() >= parameters.Q * (parameters.Q <= 1.0 ? p->getW() : 1.0)) {
627 data.push_back(*p);
628 }
629 }
630 }
631
632 if (!data.empty()) {
633
634 sort(data.begin(), data.end(), make_comparator(&hit_type::getToA, JComparison::lt())); // sort according time-of-arrival
635
636 buffer_type buffer(factoryLimit); // local buffer of hits
637 event_type out[2]; // FIFO of events
638
639 TH1D* h1 = M1[i->first];
640 TH1D* h2 = M2[i->first];
641 TH1D* h3 = M3[i->first];
642 TH1D* h4 = M4[i->first];
643 TH1D* h5 = M5[i->first];
644
645 for (buffer_type::const_iterator p = data.begin(); p != data.end(); ++p) {
646
647 buffer_type::const_iterator q = p;
648
649 while (++q != data.end() && q->getToA() - p->getToA() <= parameters.TMax_s) {}
650
651 h1->Fill(distance(p,q));
652
653 if (distance(p,q) >= parameters.numberOfHits) {
654
655 if (distance(p,q) < factoryLimit) {
656
657 if ((int) out[1].size() >= parameters.numberOfHits && velo(out[1].getVertex(), p, q) >= parameters.numberOfHits) {
658
659 // ?
660
661 } else {
662
663 buffer_type::iterator root = buffer.begin();
664 buffer_type::iterator __p = buffer.begin();
665 buffer_type::iterator __q = buffer.begin();
666
667 *root = *p;
668
669 ++__p;
670 ++__q;
671
672 for (buffer_type::const_iterator i = p; ++i != q; ) {
673 if (match(*p,*i)) {
674 *__q = *i;
675 ++__q;
676 }
677 }
678
679 h2->Fill(distance(root,__q));
680
681 if (distance(root,__q) >= parameters.numberOfHits) {
682
683 __q = clusterize(__p, __q, match);
684
685 h3->Fill(distance(root,__q));
686
687 if (distance(root,__q) >= parameters.numberOfHits) {
688
689 vertex_type vx = Velo(center, *root, __p, __q, parameters.numberOfHits);
690
691 if (vx.N >= parameters.numberOfHits) {
692
693 vx = velo(vx.getPosition(), *root, __p, __q, parameters.numberOfHits);
694
695 h4->Fill(vx.N);
696
697 if (vx.N >= parameters.numberOfHits) {
698 out[1] = event_type(vx.getVertex(), JEvent(detector.getID(), out[0].getCounter() + 1, i->first, root, __q));
699 }
700 }
701 }
702 }
703 }
704
705 } else {
706
707 out[1] = event_type(JEvent(detector.getID(), out[0].getCounter() + 1, i->first, p, q));
708 }
709
710 if (!out[1].empty()) {
711
712 if (out[0].empty()) {
713
714 out[0] = out[1]; // shift
715
716 } else if (overlap2D(out[0],out[1])) {
717
718 out[0].merge(out[1]); // merge
719
720 } else {
721
722 STATUS("trigger: " << out[0] << endl);
723
724 h5->Fill(out[0].size());
725
726 queue.push_back(out[0]); // store
727
728 out[0] = out[1]; // shift
729 }
730
731 out[1].clear();
732 }
733 }
734 }
735
736 if (!out[0].empty()) {
737
738 queue.push_back(out[0]); // store
739 }
740 STATUS(endl);
741 }
742 }
743
744 if (!queue.empty()) {
745
746 sort(queue.begin(), queue.end());
747
748 for (vector<event_type>::iterator p = queue.begin(); p != queue.end(); ) {
749
750 // select event with highest quality within given time window
751
752 vector<event_type>::iterator q = p;
753 vector<event_type>::iterator i = p;
754
755 STATUS(endl);
756
757 for ( ; q != queue.end() && overlap1D(*p,*q); ++q) {
758
759 STATUS("| " << *q << endl);
760
761 if (q->getQ() > i->getQ()) {
762 i = q;
763 }
764 }
765
766 STATUS("event: " << *i << endl);
767
768 outputFile.put(*i);
769
770 G2[i->getID()].put(i->getX(), i->getY());
771
772 p = q;
773 }
774 }
775
776 JMultipleFileScanner<JMeta> io(inputFile);
777
778 io >> outputFile;
779
780 for (const auto M : { &M1, &M2, &M3, &M4, &M5 }){
781 for (const auto& i : *M) {
782 outputFile.put(*i.second);
783 }
784 }
785
786 for (map<int, JGraph_t>::const_iterator i = G2.begin(); i != G2.end(); ++i) {
787 outputFile.put(JGraph(i->second, MAKE_CSTRING("G[" << FILL(2,'0') << i->first << "]")));
788 }
789
790 outputFile.close();
791}
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
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
JDETECTOR::JLocation::getString
int getString() const
Get string number.
Definition JLocation.hh:135
JEEP::JProperties
Utility class to parse parameter values.
Definition JProperties.hh:503
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
JGEOMETRY3D::JVector3D
Data structure for vector in three dimensions.
Definition JVector3D.hh:36
JGEOMETRY3D::JVertex3D::getVertex
const JVertex3D & getVertex() const
Get vertex.
Definition JVertex3D.hh:75
JLANG::JException
General exception.
Definition JException.hh:24
JLANG::JObjectID::getID
int getID() const
Get identifier.
Definition JObjectID.hh:50
JPARSER::JParser
Utility class to parse command line options.
Definition JParser.hh:1698
JROOT::JManager
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys.
Definition JManager.hh:47
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::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
JAANET::getVertex
JVertex3D getVertex(const Trk &track)
Get vertex.
Definition JAAnetToolkit.hh:254
JAANET::getPosition
JPosition3D getPosition(const Vec &pos)
Get position.
Definition JAAnetToolkit.hh:182
JACOUSTICS::getWaveformID
int getWaveformID(int id)
Get waveform identifier.
Definition JAcousticsToolkit.hh:39
JACOUSTICS::hydrophones_container
JContainer< std::vector< JHydrophone > > hydrophones_container
Definition JSydney.cc:80
JACOUSTICS::getSoundVelocity
static const JSoundVelocity getSoundVelocity(1541.0, -17.0e-3, -2000.0)
Function object for velocity of sound.
JDETECTOR::load
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Definition JDetectorToolkit.hh:536
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
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
FILL
Auxiliary data structure for sequence of same character.
Definition JManip.hh:330
JAANET::detector
Detector file.
Definition JHead.hh:227
JACOUSTICS::JAbstractSoundVelocity
Interface for depth dependend velocity of sound.
Definition JAbstractSoundVelocity.hh:20
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::compare
Auxiliary class to compare transmissions.
Definition JTransmission.hh:230
JACOUSTICS::JTransmission::equals
Auxiliary class to compare transmissions.
Definition JTransmission.hh:198
JACOUSTICS::JTransmission
Acoustic transmission.
Definition JTransmission.hh:31
JACOUSTICS::JTriggerParameters::Q
double Q
minimal quality if larger than one; else minimal normalised quality
Definition JAcoustics/JTriggerParameters.hh:47
JACOUSTICS::JTriggerParameters::TMax_s
double TMax_s
maximal difference between times of emission [s]
Definition JAcoustics/JTriggerParameters.hh:48
JACOUSTICS::JTriggerParameters::numberOfHits
int numberOfHits
minimal number of hits to trigger event
Definition JAcoustics/JTriggerParameters.hh:49
JACOUSTICS::event_type
Event with vertex.
Definition JAcousticsTriggerProcessor.cc:143
JACOUSTICS::event_type::TOA_s
static double TOA_s
start time of data
Definition JAcousticsTriggerProcessor.cc:217
JACOUSTICS::hit_type
Transmission with position.
Definition JBillabong.cc:70
JACOUSTICS::vertex_type
Vertex.
Definition JAcousticsTriggerProcessor.cc:266
JACOUSTICS::vertex_type::N
int N
number of hits
Definition JAcousticsTriggerProcessor.cc:292
JEEP::JContainer
Auxiliary wrapper for I/O of container with optional comment (see JComment).
Definition JContainer.hh:42
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
JROOT::JGraph
Auxiliary data structure to build TGraph.
Definition JGraph.hh:44
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::JHashMap
General purpose class for hash map of unique keys.
Definition JHashMap.hh:75
Generated by doxygen 1.12.0