JAcousticsTriggerProcessor.cc File Reference

Main program to trigger acoustic data. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include <map>
#include <algorithm>
#include "TROOT.h"
#include "TFile.h"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JHydrophone.hh"
#include "JLang/JPredicate.hh"
#include "JLang/JComparator.hh"
#include "JTools/JRange.hh"
#include "JTools/JQuantile.hh"
#include "JTools/JHashMap.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMeta.hh"
#include "JAcoustics/JToA.hh"
#include "JAcoustics/JEmitterID.hh"
#include "JAcoustics/JTransmission.hh"
#include "JAcoustics/JReceiver.hh"
#include "JAcoustics/JSoundVelocity.hh"
#include "JAcoustics/JAcousticsToolkit.hh"
#include "JAcoustics/JAcousticsSupportkit.hh"
#include "JAcoustics/JTriggerParameters.hh"
#include "JAcoustics/JEvent.hh"
#include "JAcoustics/JSupport.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
	Acoustic hit. More...
class	JACOUSTICS::JMatch3D
	3D match criterion for acoustic signals. More...
class	JACOUSTICS::JEventOverlap
	Match of two events considering overlap in time. More...

Namespaces
	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 182 of file JAcousticsTriggerProcessor.cc.

{
  using namespace std;
  using namespace JPP;
  
  typedef JTYPELIST<JEvent, JTriggerParameters, JMeta>::typelist  typelist;
 
  typedef JContainer< vector<JHydrophone> >  hydrophones_container;
 
  JMultipleFileScanner<JToA> inputFile;
  JLimit_t&                numberOfEvents = inputFile.getLimit();
  JTriggerParameters       parameters;
  int                      factoryLimit   = 10000;
  double                   TMaxExtra_s    = 100.0e-6;
  JFileRecorder<typelist>  outputFile;
  JSoundVelocity           V = getSoundVelocity;   // default sound velocity
  string                   detectorFile;
  hydrophones_container    hydrophones;            // hydrophones
  double                   precision;
  int                      debug;
 
  try {
 
    JProperties properties;
 
    properties.insert(gmake_property(parameters.Q));
    properties.insert(gmake_property(parameters.TMax_s));
    properties.insert(gmake_property(parameters.numberOfHits));
    properties.insert(gmake_property(factoryLimit));
    properties.insert(gmake_property(TMaxExtra_s));
    
    JParser<> zap("Main program to trigger acoustic data.");
    
    zap['f'] = make_field(inputFile,    "output of JConvertDB -q toashort");
    zap['n'] = make_field(numberOfEvents)                                    = JLimit::max();    
    zap['@'] = make_field(properties,   "trigger parameters");
    zap['o'] = make_field(outputFile,   "output file")                       = "event.root";
    zap['V'] = make_field(V,            "sound velocity")                    = JPARSER::initialised();
    zap['a'] = make_field(detectorFile, "detector file");
    zap['H'] = make_field(hydrophones,  "hydrophone data")                   = JPARSER::initialised();
    zap['p'] = make_field(precision,    "precision time-of-arrival")         = 1.0e-6;
    zap['W'] = make_field(getEmitterID, "waveform identification data")      = JPARSER::initialised();
    zap['d'] = make_field(debug)                                             = 1;
 
    zap(argc, argv);
  }
  catch(const exception &error) {
    FATAL(error.what() << endl);
  }
 
 
  JDetector detector;
 
  try {
    load(detectorFile, detector);
  }
  catch(const JException& error) {
    FATAL(error);
  }
 
  V.set(detector.getUTMZ());
 
 
 
  const JMatch3D      match(V, TMaxExtra_s);
  const JEventOverlap overlap(parameters.TMax_s);
 
    
  JHashMap<int, JReceiver>   receivers;
 
  for (JDetector::const_iterator i = detector.begin(); i != detector.end(); ++i) {
 
    JPosition3D pos(0.0, 0.0, 0.0);
 
    if (i->getFloor() == 0) {                    // get relative position of hydrophone
 
      try {
        pos = getPosition(hydrophones.begin(),
                          hydrophones.end(),
                          make_predicate(&JHydrophone::getString, i->getString()));
      }
      catch(const exception&) {
        continue;                                // if no position available, discard hydrophone
      }
    }
 
    receivers[i->getID()] = JReceiver(i->getID(), 
                                      i->getPosition() + pos, 
                                      i->getT0() * 1.0e-9);
  }
 
  outputFile.open();
 
  if (!outputFile.is_open()) {
    FATAL("Error opening file " << outputFile << endl);
  }
 
  outputFile.put(JMeta(argc, argv));
  outputFile.put(parameters);
 
 
  // input data
 
  typedef vector<hit_type> buffer_type;          // data type
 
  map<int, map< int, buffer_type > >  f1;        // emitter -> receiver -> data
 
  while (inputFile.hasNext()) {
 
    if (inputFile.getCounter()%100000 == 0) {
      STATUS("counter: " << setw(8) << inputFile.getCounter() << '\r' << flush); DEBUG(endl);
    }
 
    JToA* parameters = inputFile.next();
 
    if (detector.getID() != parameters->DETID) { // consistency check
      FATAL("Invalid detector identifier " << parameters->DETID << " != " << detector.getID() << endl);
    }
 
    if (receivers.has(parameters->DOMID)) {
 
      const JReceiver&    receiver = receivers[parameters->DOMID];
 
      double toa = parameters->TOA_S();
 
      const JTransmission transmission(parameters->RUN,
                                       parameters->DOMID,
                                       parameters->QUALITYFACTOR,
                                       parameters->QUALITYNORMALISATION,
                                       receiver.getT(toa),
                                       receiver.getT(toa));
      try {
        f1[getEmitterID(parameters->WAVEFORMID)][receiver.getID()].push_back(hit_type(receiver.getPosition(), transmission));
      }
      catch(const exception&) {}    
    }
  }
  STATUS(endl);
 
  for (map<int, map< int, buffer_type> >::iterator i = f1.begin(); i != f1.end(); ++i) {
 
    buffer_type data;
 
    for (map< int, buffer_type>::iterator receiver = i->second.begin(); receiver != i->second.end(); ++receiver) {
 
      // filter similar hits
 
      sort(receiver->second.begin(), receiver->second.end(), JTransmission::compare(precision));
 
      buffer_type::iterator __end = unique(receiver->second.begin(), receiver->second.end(), JTransmission::equals(precision));
 
      // selection based on quality
 
      for (buffer_type::const_iterator p = receiver->second.begin(); p != __end; ++p) {
        if (p->getQ() >= parameters.Q * (parameters.Q <= 1.0 ? p->getW() : 1.0)) {
          data.push_back(*p);
        }
      }
    }
 
    sort(data.begin(), data.end(), make_comparator(&hit_type::getToA, JComparison::lt()));          // sort according time-of-arrival
 
    buffer_type buffer(factoryLimit);               // local buffer of hits
    JEvent      out[2];                              // FIFO of events
 
    for (buffer_type::const_iterator p = data.begin(); p != data.end(); ++p) {
 
      if (distance(data.cbegin(),p)%1000 == 0) {
        STATUS("processed[" << i->first << "]: " << FIXED(5,1) << (double) (100 * distance(data.cbegin(),p)) / (double) data.size() << "%" << '\r' << flush); DEBUG(endl);
      }
 
      buffer_type::const_iterator q = p;
 
      while (++q != data.end() && q->getToA() - p->getToA() <= parameters.TMax_s) {}
 
      if (distance(p,q) >= parameters.numberOfHits) {
 
        if (distance(p,q) < factoryLimit) {
 
          buffer_type::iterator root = buffer.begin();
          buffer_type::iterator __p  = buffer.begin();
          buffer_type::iterator __q  = buffer.begin();
 
          *root = *p;
        
          ++__p;
          ++__q;
 
          for (buffer_type::const_iterator i = p; ++i != q; ) {
            if (match(*p,*i)) {
              *__q = *i;
              ++__q;
            }
          }
 
          if (distance(root,__q) >= parameters.numberOfHits) { 
        
            __q = clusterize(__p, __q, match, parameters.numberOfHits - 1);
 
            if (distance(root,__q) >= parameters.numberOfHits) {
              out[1] = JEvent(detector.getID(), out[0].getCounter() + 1, i->first, root, __q);
            }
          }
          
        } else {
 
          out[1] = JEvent(detector.getID(), out[0].getCounter() + 1, i->first, p, q);
        }
 
        if (out[0].empty()) {
 
          out[0] = out[1];                           // shift
 
        } else if (overlap(out[0],out[1])) {
      
          out[0].merge(out[1]);                      // merge
 
        } else {
            
          outputFile.put(out[0]);                    // write
 
          out[0] = out[1];                           // shift
        }
 
        out[1].clear();
      }
    }
 
    if (!out[0].empty()) {
 
      outputFile.put(out[0]);                        // write
    }
    STATUS(endl);
    STATUS("triggers[" << i->first << "]: " << setw(7) << out[0].getCounter() << endl);
  }
 
  JMultipleFileScanner<JMeta> io(inputFile);
 
  io >> outputFile;
 
  outputFile.close();
}