Auxiliary program to generate noise in Monte Carlo event. More...

#include <string>
#include <iostream>
#include <iomanip>
#include <vector>
#include "TROOT.h"
#include "TFile.h"
#include "TRandom3.h"
#include "km3net-dataformat/offline/Head.hh"
#include "km3net-dataformat/offline/MultiHead.hh"
#include "km3net-dataformat/offline/Evt.hh"
#include "JAAnet/JAAnetToolkit.hh"
#include "JSirene/JSirene.hh"
#include "JDetector/JDetector.hh"
#include "JDetector/JDetectorToolkit.hh"
#include "JDetector/JK40DefaultSimulator.hh"
#include "JDetector/JTimeRange.hh"
#include "JSupport/JMultipleFileScanner.hh"
#include "JSupport/JFileRecorder.hh"
#include "JSupport/JMonteCarloFileSupportkit.hh"
#include "JSupport/JSupport.hh"
#include "JSupport/JMeta.hh"
#include "Jeep/JParser.hh"
#include "Jeep/JMessage.hh"

Functions
int	main (int argc, char **argv)

Detailed Description

Auxiliary program to generate noise in Monte Carlo event.

This program will be deprecated (JTriggerEfficiency.cc should be used instead).

Author: mdejong

Definition in file JModK40.cc.

Function Documentation

int main	(	int	argc,
		char **	argv
	)

Definition at line 40 of file JModK40.cc.

 {
   using namespace std;
   using namespace JPP;
 
   JMultipleFileScanner<Evt>                                              inputFile;
   JFileRecorder<JTYPELIST<JAAnetTypes_t, JMeta, JRootTypes_t>::typelist> outputFile;
   JLimit_t&      numberOfEvents = inputFile.getLimit();
   string         detectorFile;
   JK40Rates      rates_Hz;
   JTimeRange     period;
   UInt_t         seed;
   int            debug;
 
   try { 
 
     JParser<> zap("Auxiliary program to generate noise in Monte Carlo event.");
     
     zap['f'] = make_field(inputFile)       = JPARSER::initialised();
     zap['o'] = make_field(outputFile)      = "modk40.root";
     zap['n'] = make_field(numberOfEvents)  = JLimit::max();
     zap['a'] = make_field(detectorFile);
     zap['T'] = make_field(period)          = JTimeRange::DEFAULT_RANGE();
     zap['B'] = make_field(rates_Hz)        = JPARSER::initialised();
     zap['S'] = make_field(seed)            = 0;
     zap['d'] = make_field(debug)           = 1;
 
     zap(argc, argv);
   }
   catch(const exception &error) {
     FATAL(error.what() << endl);
   }
 
 
   gRandom->SetSeed(seed);
 
   const int NPE = 1;
 
  
   JDetector detector;
 
   try {
     load(detectorFile, detector);
   }
   catch(const JException& error) {
     FATAL(error);
   }
 
 
   if (period == JTimeRange::DEFAULT_RANGE()) {
 
     const double Tmax = getMaximalTime(detector);
 
     period = JTimeRange(-Tmax, +Tmax);
   }
 
 
   const JK40DefaultSimulator modk40(rates_Hz);
 
   outputFile.open();
 
   if (!outputFile.is_open()) {
     FATAL("Error opening file " << outputFile << endl);
   }
 
   outputFile.put(JMeta(argc, argv));
   outputFile.put(*gRandom);
 
   if (!inputFile.empty()) {
 
     JMultipleFileScanner<Head> io(inputFile);
 
     io >> outputFile;
 
     while (inputFile.hasNext()) {
 
       STATUS("event: " << setw(10) << inputFile.getCounter() << '\r'); DEBUG(endl);
 
       Evt* event = inputFile.next();
 
       JTimeRange time_range(JTimeRange::DEFAULT_RANGE());
 
       vector<Hit>::iterator __end = event->mc_hits.end();
 
       for (vector<Hit>::iterator i = event->mc_hits.begin(); i != __end; ) {
 
         if (is_noise(*i)) {
 
           iter_swap(i, --__end);
 
         } else {
 
           time_range.include(getTime(*i));
         
           ++i;
         }
       }
 
       event->mc_hits.erase(__end, event->mc_hits.end());
     
       if (time_range.is_valid())
         time_range.add(period);
       else
         time_range = period;
 
 
       JModuleData buffer;
 
       for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
 
         buffer.reset(module->size());
       
         modk40.generateHits(*module, time_range, buffer);
       
         for (unsigned int pmt = 0; pmt != buffer.size(); ++pmt) {
         
           const JModuleData::value_type& frame = buffer[pmt];
         
           for (JModuleData::value_type::const_iterator hit = frame.begin(); hit != frame.end(); ++hit) {
           
             event->mc_hits.push_back(JHit_t(event->mc_hits.size() + 1,
                                             module->getPMT(pmt).getID(),
                                             HIT_TYPE_NOISE,
                                             0,
                                             hit->t_ns,
                                             NPE));
           }
         }
       }
 
 
       outputFile.put(*event);
     }
   } else {
 
     for (counter_type counter = 0; counter != numberOfEvents; ++counter) {
 
       STATUS("event: " << setw(10) << counter << '\r'); DEBUG(endl);
 
       Evt event;
 
       JModuleData buffer;
 
       for (JDetector::const_iterator module = detector.begin(); module != detector.end(); ++module) {
 
         buffer.reset(module->size());
       
         modk40.generateHits(*module, period, buffer);
       
         for (unsigned int pmt = 0; pmt != buffer.size(); ++pmt) {
         
           const JModuleData::value_type& frame = buffer[pmt];
         
           for (JModuleData::value_type::const_iterator hit = frame.begin(); hit != frame.end(); ++hit) {
           
             event.mc_hits.push_back(JHit_t(event.mc_hits.size() + 1,
                                            module->getPMT(pmt).getID(),
                                            HIT_TYPE_NOISE,
                                            0,
                                            hit->t_ns,
                                            NPE));
           }
         }
       }
 
 
       outputFile.put(event);
     }
   }
   STATUS(endl);
 
   outputFile.put(*gRandom);
   outputFile.close();
 }

Functions

Detailed Description

Function Documentation