Jpp  18.0.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Functions
JFitK40.cc File Reference

Auxiliary program to fit PMT parameters from JMergeCalibrateK40.cc output. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <map>
#include <cmath>
#include "km3net-dataformat/online/JDAQ.hh"
#include "km3net-dataformat/online/JDAQHeader.hh"
#include "TROOT.h"
#include "TFile.h"
#include "TH1D.h"
#include "TH2D.h"
#include "TF2.h"
#include "TMath.h"
#include "TFitResult.h"
#include "Math/MinimizerOptions.h"
#include "JPhysics/JConstants.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JPMTParametersMap.hh"
#include "JDetector/JPMTParametersToolkit.hh"
#include "JROOT/JRootFileWriter.hh"
#include "JROOT/JRootToolkit.hh"
#include "JTools/JRange.hh"
#include "JTools/JQuantile.hh"
#include "JSupport/JMeta.hh"
#include "JSupport/JSingleFileScanner.hh"
#include "JCalibrate/JCalibrateK40.hh"
#include "JCalibrate/JFitK40.hh"
#include "JCalibrate/JTDC_t.hh"
#include "Jeep/JProperties.hh"
#include "Jeep/JPrint.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Go to the source code of this file.

Functions

int main (int argc, char **argv)
 

Detailed Description

Auxiliary program to fit PMT parameters from JMergeCalibrateK40.cc output.


To force drawing of fit functions bin-by-bin compatible to the corresponding histogram data,
the function values will be stored as histograms, rather than as associated function objects.
See JCalibrateK40.hh for the corresponding extensions of histogram names.

Author
mdejong

Definition in file JFitK40.cc.

Function Documentation

int main ( int  argc,
char **  argv 
)

Definition at line 77 of file JFitK40.cc.

78 {
79  using namespace std;
80  using namespace JPP;
81  using namespace KM3NETDAQ;
82 
83  typedef JRange<double> JRange_t;
84 
85  const double epsilon = 1.0e-10; // minimal error for ROOT plotting
86  JFitK40Parameters fitk40; // setting of internal fit parameters
87 
88  string inputFile;
89  string outputFile;
90  string detectorFile;
91  string pmtFile;
92  JTDC_t TDC;
93  bool reverse;
94  bool overwriteDetector;
95  bool writeFits;
96  bool fitAngDist;
97  bool fitModel;
98  JRange_t X;
99  JRange_t Y;
100  string option;
101  int debug;
102 
103  try {
104 
105  JProperties properties;
106 
107  properties.insert(gmake_property(STDEV));
108  properties.insert(gmake_property(QE_MIN_VALUE));
109  properties.insert(gmake_property(QE_MAX_ERROR));
110  properties.insert(gmake_property(MINIMAL_RATE_HZ));
111  properties.insert(gmake_property(fitk40.Rate_Hz));
112  properties.insert(gmake_property(fitk40.p1));
113  properties.insert(gmake_property(fitk40.p2));
114  properties.insert(gmake_property(fitk40.p3));
115  properties.insert(gmake_property(fitk40.p4));
116 
117  JParser<> zap("Auxiliary program to fit PMT parameters from JMergeCalibrateK40 output.");
118 
119  zap['@'] = make_field(properties) = JPARSER::initialised();
120  zap['f'] = make_field(inputFile, "input file (output from JMergeCalibrateK40).");
121  zap['o'] = make_field(outputFile, "output file.") = "fit.root";
122  zap['a'] = make_field(detectorFile, "detector file.");
123  zap['P'] = make_field(pmtFile, "specify PMT file name that can be given to JTriggerEfficiency.") = "";
124  zap['!'] = make_field(TDC,
125  "Fix time offset(s) of PMT(s) of certain module(s), e.g."
126  "\n-! \"808969848 0 808982077 23\" will fix offset of PMT 0 of module 808969848 and of PMT 23 of module 808982077."
127  "\nSame PMT offset can be fixed for all optical modules, e.g."
128  "\n-! \"-1 0 -1 23\" will fix offsets of PMTs 0 and 23 for all optical modules.")
129  = JPARSER::initialised();
130  zap['r'] = make_field(reverse, "reverse TDC constraints due to option -! <TDC>.");
131  zap['A'] = make_field(overwriteDetector, "overwrite detector file provided through '-a' with correct time offsets.");
132  zap['w'] = make_field(writeFits, "write fit results.");
133  zap['D'] = make_field(fitAngDist, "fit angular distribution; fix normalisation.");
134  zap['M'] = make_field(fitModel, "fit angular distribution; fix PMT QEs = 1.0.");
135  zap['x'] = make_field(X, "ROOT fit range (PMT pairs).") = JRange_t();
136  zap['y'] = make_field(Y, "ROOT fit range (time residual).") = JRange_t();
137  zap['O'] = make_field(option, "ROOT fit option, see TH1::Fit.") = "";
138  zap['d'] = make_field(debug, "debug.") = 1;
139 
140  zap(argc, argv);
141  }
142  catch(const exception &error) {
143  FATAL(error.what() << endl);
144  }
145 
146 
147  ROOT::Math::MinimizerOptions::SetDefaultMaxFunctionCalls(1000000);
148 
149 
150  if (reverse) {
151  TDC.reverse();
152  }
153 
154  for (JTDC_t::const_iterator i = TDC.begin(); i != TDC.end(); ++i) {
155  DEBUG("PMT " << setw(10) << i->first << ' ' << setw(2) << i->second << " constrain t0." << endl);
156  }
157 
158  TDC.is_valid(true);
159 
160 
161  //----------------------------------------------------------
162  // load detector file and PMT parameters
163  //----------------------------------------------------------
164 
165  JDetector detector;
166 
167  try {
168  load(detectorFile, detector);
169  }
170  catch(const JException& error) {
171  FATAL(error);
172  }
173 
174 
175  JPMTParametersMap parameters;
176 
177  if (pmtFile != "") {
178  try {
179  parameters.load(pmtFile.c_str());
180  }
181  catch(const JException& error) {
182  //ERROR(error.what());
183  }
184  }
185 
186 
187  if (option.find('O') == string::npos) { option += '0'; }
188  if (option.find('S') == string::npos) { option += 'S'; }
189  if (option.find('Q') == string::npos && debug < JEEP::debug_t) { option += 'Q'; }
190 
191 
192  TFile* in = TFile::Open(inputFile.c_str(), "exist");
193 
194  if (in == NULL || !in->IsOpen()) {
195  FATAL("File: " << inputFile << " not opened." << endl);
196  }
197 
198 
199  TDirectory::AddDirectory(NULL);
200 
201  TFile out(outputFile.c_str(), "recreate");
202 
203  // Histograms with global fit results.
204 
205  TH1D hc("chi2", NULL, 500, 0.0, 5.0);
206 
207  TH1D p1("p1", NULL, 501, -5.0, +5.0);
208  TH1D p2("p2", NULL, 501, -5.0, +5.0);
209  TH1D p3("p3", NULL, 501, -5.0, +5.0);
210  TH1D p4("p4", NULL, 501, -5.0, +5.0);
211  TH1D bg("bg", NULL, 501, -0.1, +0.1);
212  TH1D cc("cc", NULL, 501, -0.1, +0.1);
213 
214 
215  //----------------------------------------------------------
216  // loop over modules in detector file to fit histogram
217  //----------------------------------------------------------
218 
219  for (JDetector::iterator module = detector.begin(); module != detector.end(); ++module) {
220 
221  NOTICE("Module " << setw(10) << module->getID() << ' ' << module->getLocation() << endl);
222 
223  if (module->empty()) {
224  continue;
225  }
226 
227  //----------------------------------------------------------
228  // get histogram for this module
229  //----------------------------------------------------------
230 
231  TH2D* h2s = (TH2D*) in->Get(MAKE_CSTRING(module->getID() << _2R));
232 
233  if (h2s == NULL || h2s->GetEntries() == 0) {
234 
235  NOTICE("No data for module " << module->getID() << "; set corresponding QEs to 0." << endl);
236 
237  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
238  parameters[JPMTIdentifier(module->getID(), pmt)].QE = 0.0;
239  }
240 
241  continue;
242  }
243 
244  //----------------------------------------------------------
245  // get constraints
246  //----------------------------------------------------------
247 
248  const JTDC_t::range_type range = TDC.equal_range(module->getID());
249 
250  //----------------------------------------------------------
251  // set-up fitting procedure
252  //----------------------------------------------------------
253 
254  JFitK40 fit(*module,
255  h2s->GetXaxis()->GetXmin(),
256  h2s->GetXaxis()->GetXmax(),
257  h2s->GetYaxis()->GetXmin(),
258  h2s->GetYaxis()->GetXmax(),
259  range.first == range.second);
260 
261  fit.setModelParameters(fitk40.getModelParameters());
262 
263  for (JTDC_t::const_iterator i = range.first; i != range.second; ++i) {
264  fixParameter(fit, fit.getModelParameter(i->second, &JPMTParameters_t::t0));
265  }
266 
267  if (fitModel) {
268  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
269  fixParameter(fit, fit.getModelParameter(pmt, &JPMTParameters_t::QE));
270  }
271 
272  } else {
273 
274  fixParameter(fit, fit.getModelParameter(&JFitK40::Rate_Hz));
275 
276  if (!fitAngDist) {
277  fixParameter(fit, fit.getModelParameter(&JFitK40::p1));
278  fixParameter(fit, fit.getModelParameter(&JFitK40::p2));
279  fixParameter(fit, fit.getModelParameter(&JFitK40::p3));
280  fixParameter(fit, fit.getModelParameter(&JFitK40::p4));
281  }
282  }
283 
284  //----------------------------------------------------------
285  // check validity of data
286  //----------------------------------------------------------
287 
288  vector<int> buffer(NUMBER_OF_PMTS, 1);
289 
290  for (int ix = 1; ix <= h2s->GetXaxis()->GetNbins(); ++ix) {
291 
292  Double_t value = 0.0;
293  Double_t error = 0.0;
294 
295  for (int iy = 1; iy <= h2s->GetYaxis()->GetNbins(); ++iy) {
296  value += h2s->GetBinContent(ix,iy);
297  error += h2s->GetBinError(ix,iy) * h2s->GetBinError(ix,iy);
298  }
299 
300  error = sqrt(error);
301 
302  if (value < MINIMAL_RATE_HZ + STDEV*error) {
303 
304  const JFitK40::pair_type pair = fit.getPair(ix-1);
305 
306  buffer[pair.first] += 1;
307  buffer[pair.second] += 1;
308  }
309  }
310 
311  try {
312 
313  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
314 
315  if (buffer[pmt] == NUMBER_OF_PMTS) {
316 
317  WARNING("PMT " << setw(10) << module->getID() << ' ' << setw(2) << pmt << " too low a rate; switch off." << endl);
318 
319  fit.disablePMT(pmt);
320  }
321  }
322  }
323  catch(const JException& error) {
324 
325  NOTICE("Skip fit for module " << module->getID() << "; set corresponding QEs to 0." << endl);
326 
327  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
328  parameters[JPMTIdentifier(module->getID(), pmt)].QE = 0.0;
329  }
330 
331  continue;
332  }
333 
334  // constrain fit to specified range(s)
335 
336  if (X != JRange_t()) {
337  h2s->GetXaxis()->SetRangeUser(max(X.getLowerLimit(), h2s->GetXaxis()->GetXmin()),
338  min(X.getUpperLimit(), h2s->GetXaxis()->GetXmax()));
339  }
340  if (Y != JRange_t()) {
341  h2s->GetYaxis()->SetRangeUser(max(Y.getLowerLimit(), h2s->GetYaxis()->GetXmin()),
342  min(Y.getUpperLimit(), h2s->GetYaxis()->GetXmax()));
343  }
344 
345  // 1D-histogram with peak positions
346 
347  TH1D h1(MAKE_CSTRING(module->getID() << _1D), NULL, h2s->GetXaxis()->GetNbins(), h2s->GetXaxis()->GetXmin(), h2s->GetXaxis()->GetXmax());
348 
349  const Double_t sigma[] = {
350 
351  (min(Y.getUpperLimit(), h2s->GetYaxis()->GetXmax()) -
352  max(Y.getLowerLimit(), h2s->GetYaxis()->GetXmin())) / sqrt(12.0),
353 
354  sqrt(JPMTParameters_t().TTS * JPMTParameters_t().TTS * 2.0 + fit.getSigmaK40() * fit.getSigmaK40())
355  };
356 
357 
358  for (Int_t ix = 1; ix <= h2s->GetXaxis()->GetNbins(); ++ix) {
359 
360  Double_t y0 = 0.0;
361  Double_t zmax = 0.0;
362  Double_t value = 0.0;
363  Double_t error = 0.0;
364 
365  for (Int_t iy = 1; iy <= h2s->GetYaxis()->GetNbins(); ++iy) {
366 
367  const Double_t y = h2s->GetYaxis()->GetBinCenter(iy);
368  const Double_t z = h2s->GetBinContent(ix, iy);
369  const Double_t w = h2s->GetBinError (ix, iy);
370 
371  if (Y(y)) {
372 
373  if (z > zmax) {
374  y0 = y;
375  zmax = z;
376  }
377 
378  value += z;
379  error += w*w;
380  }
381  }
382 
383  error = sqrt(error);
384 
385  h1.SetBinContent(ix, y0);
386  h1.SetBinError (ix, value > MINIMAL_RATE_HZ + STDEV * error ? sigma[1] * error / value : sigma[0]);
387  }
388 
389  {
390  JFitK40_t<TF1> f1(fit);
391 
392  TFitResultPtr result = f1(h1, option.c_str());
393 
394  if (result.Get() == NULL) {
395  FATAL("Invalid TFitResultPtr " << h1.GetName() << endl);
396  }
397 
398  if (writeFits) {
399 
400  TH1D h2(MAKE_CSTRING(module->getID() << _1F), NULL, h1.GetXaxis()->GetNbins(), h1.GetXaxis()->GetXmin(), h1.GetXaxis()->GetXmax());
401 
402  for (int ix = 1; ix <= h2.GetXaxis()->GetNbins(); ++ix) {
403 
404  const Double_t x = h2.GetXaxis()->GetBinCenter(ix);
405 
406  h2.SetBinContent(ix, f1.Eval(x));
407  h2.SetBinError (ix, 0.0);
408  }
409 
410  out << h1 << h2;
411  }
412 
413  if (debug >= JEEP::notice_t) {
414 
415  cout << "Fit result chi2 / NDF " << result->Chi2() << ' ' << result->Ndf() << ' ' << (result->IsValid() ? "" : "failed") << endl;
416 
417  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
418 
419  cout << ' '
420  << setw(2) << pmt << ' '
421  << FIXED(7,2) << f1.getT0 (pmt) << ' ';
422  cout << (result->IsParameterFixed(f1.getModelParameter(pmt, &JPMTParameters_t::t0)) ? " *** fixed t0 *** " : "");
423  cout << (f1.is_disabled (pmt) ? " (disabled)" : "");
424  cout << (f1.is_average_t0(pmt) ? " (averaged)" : "");
425  cout << endl;
426  }
427  }
428 
429  fit.setModelParameters(f1.getModelParameters());
430  }
431 
432 
433  TFitResultPtr result = fit(*h2s, option);
434 
435 
436  bool refit = false;
437 
438  try {
439 
440  const JFitK40Parameters values(fit.GetParameters());
441  const JFitK40Parameters errors(fit.GetParErrors());
442 
443  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
444 
445  if (!result->IsParameterFixed(fit.getModelParameter(pmt, &JPMTParameters_t::QE)) && !is_valid(values.getQE(pmt), errors.getQE(pmt))) {
446 
447  WARNING("PMT " << setw(10) << module->getID() << ' ' << setw(2) << pmt << ' '
448  << "measured QE "
449  << FIXED(5,2) << values.getQE(pmt) << " +/- "
450  << FIXED(5,2) << errors.getQE(pmt)
451  << " compatible with zero -> set QE = 0 and t0 = 0; -> refit." << endl);
452 
453 
454  fit.disablePMT(pmt);
455 
456  refit = true;
457  }
458  }
459  }
460  catch(const JException& error) {
461 
462  NOTICE("Skip fit for module " << module->getID() << "; set corresponding QEs to 0." << endl);
463 
464  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
465  parameters[JPMTIdentifier(module->getID(), pmt)].QE = 0.0;
466  }
467 
468  continue;
469  }
470 
471 
472  if (refit) {
473 
474  JFitK40_t<TF1> f1(fit);
475 
476  f1(h1, option);
477 
478  fit.setModelParameters(f1.getModelParameters());
479 
480  result = fit(*h2s, option);
481  }
482 
483 
484  const JFitK40Parameters values(fit.GetParameters());
485  const JFitK40Parameters errors(fit.GetParErrors());
486 
487  hc.Fill(TMath::Log10(result->Chi2() / result->Ndf()));
488 
489  TH2D h2f(MAKE_CSTRING(module->getID() << _2F),
490  NULL,
491  h2s->GetXaxis()->GetNbins(), h2s->GetXaxis()->GetXmin (), h2s->GetXaxis()->GetXmax (),
492  h2s->GetYaxis()->GetNbins(), h2s->GetYaxis()->GetXmin (), h2s->GetYaxis()->GetXmax ());
493 
494  for (int ix = 1; ix <= h2f.GetXaxis()->GetNbins(); ++ix) {
495  for (int iy = 1; iy <= h2f.GetYaxis()->GetNbins(); ++iy) {
496 
497  const Double_t x = h2f.GetXaxis()->GetBinCenter(ix);
498  const Double_t y = h2f.GetYaxis()->GetBinCenter(iy);
499 
500  h2f.SetBinContent(ix, iy, fit.Eval(x,y));
501  h2f.SetBinError (ix, iy, 0.0);
502  }
503  }
504 
505 
506  if (debug >= JEEP::notice_t) {
507 
508  cout << "Fit result chi2 / NDF " << result->Chi2() << ' ' << result->Ndf() << ' ' << (result->IsValid() ? "" : "failed") << endl;
509 
510  cout << "Rate_Hz= " << FIXED(8,4) << values.Rate_Hz << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::Rate_Hz)) ? " *** fixed *** " : "") << endl;
511 
512  cout << "p1= " << FIXED(8,4) << values.p1 << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::p1)) ? " *** fixed *** " : "") << endl;
513  cout << "p2= " << FIXED(8,4) << values.p2 << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::p2)) ? " *** fixed *** " : "") << endl;
514  cout << "p3= " << FIXED(8,4) << values.p3 << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::p2)) ? " *** fixed *** " : "") << endl;
515  cout << "p4= " << FIXED(8,4) << values.p4 << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::p3)) ? " *** fixed *** " : "") << endl;
516 
517  cout << "Background (correlated): " << FIXED(8,5) << values.bg << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::bg)) ? " *** fixed *** " : "") << endl;
518  cout << "Background (uncorrelated): " << FIXED(8,5) << values.cc << (result->IsParameterFixed(fit.getModelParameter(&JFitK40::cc)) ? " *** fixed *** " : "") << endl;
519 
520  cout << " PMT t0 [ns] TTS [ns] R" << endl;
521 
522  JQuantile Q[3];
523 
524  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
525 
526  cout << ' '
527  << setw(2) << pmt << ' '
528  << FIXED(7,2) << fit.getT0 (pmt) << ' '
529  << FIXED(7,2) << fit.getTTS(pmt) << ' '
530  << FIXED(7,3) << fit.getQE (pmt);
531 
532  cout << (result->IsParameterFixed(fit.getModelParameter(pmt, &JPMTParameters_t::QE)) ? " *** fixed QE *** " : "");
533  cout << (result->IsParameterFixed(fit.getModelParameter(pmt, &JPMTParameters_t::TTS)) ? " *** fixed TTS *** " : "");
534  cout << (result->IsParameterFixed(fit.getModelParameter(pmt, &JPMTParameters_t::t0)) ? " *** fixed t0 *** " : "");
535  cout << (fit.is_disabled (pmt) ? " (disabled)" : "");
536  cout << (fit.is_average_t0(pmt) ? " (averaged)" : "");
537  cout << endl;
538 
539  Q[0].put(fit.getT0 (pmt));
540  Q[1].put(fit.getTTS(pmt));
541  Q[2].put(fit.getQE (pmt));
542  }
543 
544  cout << "<> "
545  << FIXED(7,2) << Q[0].getMean() << ' '
546  << FIXED(7,2) << Q[1].getMean() << ' '
547  << FIXED(7,3) << Q[2].getMean() << endl;
548  cout << "+/- "
549  << FIXED(7,2) << Q[0].getSTDev() << ' '
550  << FIXED(7,2) << Q[1].getSTDev() << ' '
551  << FIXED(7,3) << Q[2].getSTDev() << endl;
552  }
553 
554  h2s->Write();
555  h2f .Write();
556 
557  if (writeFits) {
558 
559  // Histograms with fit results.
560 
561  p1.Fill(values.p1);
562  p2.Fill(values.p2);
563  p3.Fill(values.p3);
564  p4.Fill(values.p4);
565  bg.Fill(values.bg);
566  cc.Fill(values.cc);
567 
568  TH1D h1t(MAKE_CSTRING(module->getID() << ".1t0"), NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS - 0.5);
569  TH1D h1s(MAKE_CSTRING(module->getID() << ".1TTS"), NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS - 0.5);
570  TH1D h1q(MAKE_CSTRING(module->getID() << ".1QE"), NULL, NUMBER_OF_PMTS, -0.5, NUMBER_OF_PMTS - 0.5);
571 
572  h1t.Sumw2();
573  h1s.Sumw2();
574  h1q.Sumw2();
575 
576  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
577 
578  h1t.SetBinContent(pmt + 1, values.getT0 (pmt));
579  h1s.SetBinContent(pmt + 1, values.getTTS(pmt));
580  h1q.SetBinContent(pmt + 1, values.getQE (pmt));
581 
582  h1t.SetBinError (pmt + 1, max(errors.getT0 (pmt), epsilon));
583  h1s.SetBinError (pmt + 1, max(errors.getTTS(pmt), epsilon));
584  h1q.SetBinError (pmt + 1, max(errors.getQE (pmt), epsilon));
585  }
586 
587  out << h1t << h1s << h1q;
588  }
589 
590 
591  // save output
592 
593  for (int pmt = 0; pmt != NUMBER_OF_PMTS; ++pmt) {
594 
595  JPMTParameters& data = parameters[JPMTIdentifier(module->getID(), pmt)];
596 
597  // correction to QE
598 
599  const double P = getSurvivalProbability(data);
600 
601  if (P > 0.0)
602  data.QE = values.getQE(pmt) / P;
603  else
604  data.QE = 0.0;
605 
606  if (overwriteDetector) {
607  module->getPMT(pmt).addT0(fit.getT0(pmt));
608  }
609  }
610  }
611 
612  for (JRootFileReader<JDAQHeader> in(inputFile.c_str()); in.hasNext(); ) {
613  putObject(out, *in.next());
614  }
615 
616  if (writeFits) {
617  out << hc << p1 << p2 << p3 << p4 << bg << cc;
618  }
619 
620  out.Close();
621 
622  {
623  vector<JMeta> buffer(1, JMeta(argc, argv));
624 
625  for (JSingleFileScanner<JMeta> in(inputFile); in.hasNext(); ) {
626 
627  const JMeta* meta = in.next();
628 
629  buffer.push_back(*meta);
630  }
631 
632  for (vector<JMeta>::const_reverse_iterator i = buffer.rbegin(); i != buffer.rend(); ++i) {
633  parameters.comment.add(*i);
634  detector .comment.add(*i);
635  }
636  }
637 
638  if (overwriteDetector) {
639 
640  NOTICE("Store calibration data on file " << detectorFile << endl);
641 
642  store(detectorFile, detector);
643  }
644 
645  if (pmtFile != "") {
646  parameters.store(pmtFile.c_str());
647  }
648 
649  return 0;
650 }
JSUPPORT::JMeta
Auxiliary class for ROOT I/O of application specific meta data.
Definition: JMeta.hh:70
JPARSER::JParser
Utility class to parse command line options.
Definition: JParser.hh:1514
JLANG::JException
General exception.
Definition: JException.hh:23
JCALIBRATE::pair_type
JCombinatorics::pair_type pair_type
Definition: JCalibrateK40.hh:41
JTOOLS::w
data_type w[N+1][M+1]
Definition: JPolint.hh:778
Q
Q(UTCMax_s-UTCMin_s)-livetime_s
JEEP::debug_t
debug
Definition: JMessage.hh:29
p1
TPaveText * p1
Definition: JDrawModule3D.cc:35
JTOOLS::JQuantile::getSTDev
double getSTDev() const
Get standard deviation.
Definition: JQuantile.hh:281
gmake_property
#define gmake_property(A)
macros to convert (template) parameter to JPropertiesElement object
Definition: JProperties.hh:1291
JTOOLS::JQuantile
Auxiliary data structure for running average, standard deviation and quantiles.
Definition: JQuantile.hh:43
JSON::is_valid
bool is_valid(const json &js)
Check validity of JSon data.
Definition: JSon/JSupport.hh:246
JDETECTOR::JDetector
Detector data structure.
Definition: JDetector.hh:89
JCALIBRATE::JTDC_t
Auxiliary class for TDC constraints.
Definition: JTDC_t.hh:37
JEEP::notice_t
notice
Definition: JMessage.hh:32
JCALIBRATE::JFitK40_t< TF2 >
Template specialisation of two-fold coincidence rate due to K40 and other radioactive decays...
Definition: JFitK40.hh:586
JEEP::JProperties
Utility class to parse parameter values.
Definition: JProperties.hh:496
parameters
*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
Definition: diff-Tuna.sh:38
MAKE_CSTRING
#define MAKE_CSTRING(A)
Make C-string.
Definition: JPrint.hh:136
JPARSER::initialised
Empty structure for specification of parser element that is initialised (i.e. does not require input)...
Definition: JParser.hh:83
JCALIBRATE::JFitK40_t< TF1 >
Template specialisation of two-fold coincidence rate due to K40 and other radioactive decays...
Definition: JFitK40.hh:716
Y
then fatal Wrong number of arguments fi set_variable STRING $argv[1] set_variable DETECTORXY_TXT $WORKDIR $DETECTORXY_TXT tail read X Y CHI2 RMS printf optimum n $X $Y $CHI2 $RMS awk v Y
Definition: detector-XY:fit1d.sh:33
FIXED
Auxiliary data structure for floating point format specification.
Definition: JManip.hh:446
outputFile
string outputFile
Definition: JDAQTimesliceSelector.cc:37
JTOOLS::JQuantile::getMean
double getMean() const
Get mean value.
Definition: JQuantile.hh:252
GAUSS_HERMITE::sigma
const double sigma[]
Definition: JQuadrature.cc:74
JCALIBRATE::_1F
static const char *const _1F
Name extension for 1D time offset fit.
Definition: JCalibrateK40.hh:39
GAUSS_LEGENDRE::f1
const JPolynome f1(1.0, 2.0, 3.0)
Function.
JAANET::JRange_t
Type definition of range.
Definition: JHead.hh:41
JCALIBRATE::_2F
static const char *const _2F
Name extension for 2F rate fitted.
Definition: JCalibrateK40.hh:37
JAANET::detector
Detector file.
Definition: JHead.hh:226
make_field
#define make_field(A,...)
macro to convert parameter to JParserTemplateElement object
Definition: JParser.hh:1989
JDETECTOR::store
void store(const std::string &file_name, const JDetector &detector)
Store detector to output file.
Definition: JDetectorToolkit.hh:578
NOTICE
#define NOTICE(A)
Definition: JMessage.hh:64
JTOOLS::JQuantile::put
void put(const double x, const double w=1.0)
Put value.
Definition: JQuantile.hh:133
JCALIBRATE::_2R
static const char *const _2R
Name extension for 2D rate measured.
Definition: JCalibrateK40.hh:36
JDETECTOR::JPMTParametersMap
Auxiliary class for map of PMT parameters.
Definition: JPMTParametersMap.hh:89
JROOT::putObject
bool putObject(TDirectory &dir, const TObject &object)
Write object to ROOT directory.
Definition: JRootFileWriter.hh:41
JCALIBRATE::JPMTParameters_t
Fit parameters for single PMT.
Definition: JFitK40.hh:46
FATAL
#define FATAL(A)
Definition: JMessage.hh:67
range
z range($ZMAX-$ZMIN)< $MINIMAL_DZ." fi fi typeset -Z 4 STRING typeset -Z 2 FLOOR JPlot1D -f $
Definition: module-Z:fit.sh:118
p2
p2
Definition: module-Z:fit.sh:74
JDETECTOR::load
void load(const std::string &file_name, JDetector &detector)
Load detector from input file.
Definition: JDetectorToolkit.hh:485
JDETECTOR::getSurvivalProbability
double getSurvivalProbability(const JPMTParameters &parameters)
Get model dependent probability that a one photo-electron hit survives the simulation of the PMT assu...
Definition: JPMTParametersToolkit.hh:32
JROOT::fixParameter
bool fixParameter(TF1 &f1, const JFitParameter_t &parameter)
Fix fit parameter.
Definition: JRootToolkit.hh:536
TDC
set_variable TDC
Definition: JPrintTDC.sh:20
X
no fit printf nominal n $STRING awk v X
Definition: detector-XY:fit.sh:102
JSUPPORT::JSingleFileScanner
Object reading from a list of files.
Definition: JSingleFileScanner.hh:75
values
possible values
Definition: detector-XY:run.sh:42
JDETECTOR::JPMTParameters
Data structure for PMT parameters.
Definition: JPMTParameters.hh:29
detector
do set_variable DETECTOR_TXT $WORKDIR detector
Definition: detector-XY:fit1d.sh:38
KM3NETDAQ::NUMBER_OF_PMTS
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
Definition: JDAQ.hh:26
in
then fatal Wrong number of arguments fi set_variable DETECTOR $argv[1] set_variable INPUT_FILE $argv[2] eval JPrintDetector a $DETECTOR O IDENTIFIER eval JPrintDetector a $DETECTOR O SUMMARY JAcoustics sh $DETECTOR_ID source JAcousticsToolkit sh CHECK_EXIT_CODE typeset A EMITTERS get_tripods $WORKDIR tripod txt EMITTERS get_transmitters $WORKDIR transmitter txt EMITTERS for EMITTER in
Definition: JCanberra.sh:46
WARNING
then echo WARNING
Definition: JTuneHV.sh:91
JCALIBRATE::JFitK40Parameters
Fit parameters for two-fold coincidence rate due to K40.
Definition: JFitK40.hh:66
JDETECTOR::JPMTParameters::QE
double QE
relative quantum efficiency
Definition: JPMTParameters.hh:233
JROOT::JRootFileReader
ROOT file reader.
Definition: JRootFileReader.hh:54
p3
p3
Definition: module-Z:fit.sh:74
GAUSS_LEGENDRE::epsilon
const double epsilon
Definition: JQuadrature.cc:21
JCALIBRATE::_1D
static const char *const _1D
Name extension for 1D time offset.
Definition: JCalibrateK40.hh:38
debug
int debug
debug level
Definition: archive-put-wiki-detectors.sh:92
P
then $DIR JPlotNPE PDG P
Definition: JPlotNPE-PDG.sh:62
DEBUG
#define DEBUG(A)
Message macros.
Definition: JMessage.hh:62
Generated by   doxygen 1.8.5