48 static const int NUMBER_OF_PINGS = 11;
49 static const JWaveform EMITTER_WAVEFORM(2.5e6,
60 int main(
int argc,
char **argv)
82 disable_container disable;
83 waveform_container waveforms;
88 numberOfPings[
WILDCARD] = NUMBER_OF_PINGS;
89 waveforms [
WILDCARD] = EMITTER_WAVEFORM;
93 JParser<> zap(
"Example application to test fit of model to simulated acoustic data.");
104 zap[
'u'] =
make_field(unify,
"unify weighing of pings");
108 zap[
'B'] =
make_field(background,
"background probability") = 0.0;
114 catch(
const exception &error) {
115 FATAL(error.what() << endl);
119 gRandom->SetSeed(seed);
121 if (numberOfEvents <= 0) {
FATAL(
"Invalid number of events " << numberOfEvents << endl); }
134 for (tripods_container::const_iterator
i =
tripods.begin();
i !=
tripods.end(); ++
i) {
136 i->getUTMPosition() -
detector.getUTMPosition()));
139 for (transmitters_container::const_iterator
i = transmitters.begin();
i != transmitters.end(); ++
i) {
144 catch(
const exception&) {
149 if (
detector.empty()) {
FATAL(
"No modules in detector." << endl); }
150 if (emitters.empty()) {
FATAL(
"No emitters in system." << endl); }
155 const double D_m = -
detector.getUTMZ();
168 TH1D h0(
"cpu", NULL, 100, 0.0, 1.0e3);
169 TH1D h1(
"chi2/NDF", NULL, 100, 0.0, 5.0);
175 for (
int number_of_events = 0,
count = 0; number_of_events != numberOfEvents; ++number_of_events) {
177 STATUS(
"event: " << setw(10) << number_of_events <<
'\r');
DEBUG(endl);
181 for (JDetector::const_iterator module =
detector.begin(); module !=
detector.end(); ++module) {
186 gRandom->Uniform(-2.0e-2, +2.0e-2));
190 DEBUG(
"Model" << endl << model << endl);
194 for (JDetector::iterator module =
detector.begin(); module !=
detector.end(); ++module) {
195 if (module->getFloor() != 0 && geometry.
hasLocation(module->getLocation())) {
196 module->set(geometry[module->getString()].getPosition(model.
string[module->getString()], module->getFloor()));
203 int minimum_number_of_pings = numeric_limits<int>::max();
206 minimum_number_of_pings = min(minimum_number_of_pings,
getValue(numberOfPings, emitter->getID(), numberOfPings[
WILDCARD]));
213 const int number_of_pings =
getValue(numberOfPings, emitter->getID(), numberOfPings[
WILDCARD]);
215 const double weight = (unify ? (double) minimum_number_of_pings / (
double) number_of_pings : 1.0);
217 for (
int ping_counter = 0; ping_counter != number_of_pings; ++ping_counter) {
221 const double toe_s = emitter->getID() * 1000.0 + ping_counter * 10.0 + gRandom->Uniform(-1.0, +1.0);
225 for (JDetector::const_iterator module =
detector.begin(); module !=
detector.end(); ++module) {
227 if (disable.count(
JTransmission_t(emitter->getID(), module->getID())) == 0) {
233 const double d_m =
getDistance(module->getPosition(), emitter->getPosition());
234 const double toa_s = toe_s +
V.getTime(d_m, emitter->getZ(), module->getZ());
235 const double Q = waveform.
getQ(D_m, d_m);
241 if (gRandom->Rndm() >= background)
242 t1_s = gRandom->Gaus(toa_s,
parameters.sigma_s);
247 const JHit hit(*emitter,
249 module->getLocation(),
254 DEBUG(
"hit: " << hit <<
' ' <<
FIXED(7,1) << Q << endl);
264 const chrono::steady_clock::time_point t0 = chrono::steady_clock::now();
266 const auto result = katoomba(data.begin(), data.end());
268 const chrono::steady_clock::time_point t1 = chrono::steady_clock::now();
272 cout <<
"result:" <<
' '
276 for (data_type::const_iterator hit =
result.begin; hit !=
result.end; ++hit) {
281 h0.Fill(chrono::duration_cast<chrono::milliseconds>(t1 - t0).
count());
286 const double tx = (
i->second.tx -
result.value.string [
i->first].tx) * 1.0e3;
287 const double ty = (
i->second.ty -
result.value.string [
i->first].ty) * 1.0e3;
290 H2[
i->first]->Fill(tx, ty);
295 const double t1 =
i->second.t1 -
result.value.emission[
i->first].t1;
298 H1[
i->first.getID()]->Fill(t1);
Fit function of acoustic model.
Utility class to parse command line options.
Q(UTCMax_s-UTCMin_s)-livetime_s
double getValue(const JScale_t scale)
Get numerical value corresponding to scale.
int main(int argc, char *argv[])
std::vector< event_type > data_type
JContainer< std::vector< JTransmitter > > transmitters_container
static JDetectorMechanics getMechanics
Function object to get string mechanics.
then usage else fatal Wrong number of arguments fi JCookie sh eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O CAN eval JPrintDetector a $DETECTOR O SUMMARY source JAcousticsToolkit sh expand_array RUNS for KEY in sound_velocity waveform
then usage $script< input file >[option] nPossible options count
*fatal Wrong number of arguments esac JCookie sh typeset Z DETECTOR typeset Z SOURCE_RUN typeset Z TARGET_RUN set_variable PARAMETERS_FILE $WORKDIR parameters
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Dynamic ROOT object management.
Auxiliary data structure for floating point format specification.
V(JDAQEvent-JTriggerReprocessor)*1.0/(JDAQEvent+1.0e-10)
Data structure for detector geometry and calibration.
Global fit of prameterised detector geometry to acoustics data.
double getDistance(const JFirst_t &first, const JSecond_t &second)
Get distance between objects.
Data structure for hydrophone.
static const char WILDCARD
Model for fit to acoustics data.
static const JSoundVelocity getSoundVelocity(1541.0,-17.0e-3,-2000.0)
Function object for velocity of sound.
Auxiliary class to manage set of compatible ROOT objects (e.g. histograms) using unique keys...
JContainer< std::vector< JHydrophone > > hydrophones_container
Data structure for transmitter.
Auxiliary wrapper for I/O of container with optional comment (see JComment).
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
JPosition3D getPosition(const Vec &pos)
Get position.
JACOUSTICS::JModel::emission_type emission
void Write(TDirectory &out, const bool wm=false)
Write objects to file.
JContainer< std::vector< JTripod > > tripods_container
General purpose messaging.
Implementation for depth dependend velocity of sound.
JACOUSTICS::JModel::string_type string
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Utility class to parse command line options.
Acoustic transmission identifier.
Data structure for measured coincidence rates of all pairs of PMTs in optical module.
bool hasLocation(const JLocation &location) const
Check if this detector has given location.
do set_variable DETECTOR_TXT $WORKDIR detector
Data structure for tripod.
Acoustic transmission identifier.
JKatoomba< JAbstractMinimiser > evaluator
Template definition of fit function of acoustic model.
#define DEBUG(A)
Message macros.
Data structure for optical module.