Jpp  17.2.0
the software that should make you happy
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
JRandomTimeslice.hh
Go to the documentation of this file.
1 #ifndef __JTIMESLICE__JRANDOMTIMESLICE__
2 #define __JTIMESLICE__JRANDOMTIMESLICE__
3 
4 #include <vector>
5 
6 #include <TRandom3.h>
7 
10 #include "JDetector/JPMTAddress.hh"
11 #include "JDetector/JTimeRange.hh"
15 
16 
17 /**
18  * \author mdejong
19  */
20 
21 namespace KM3NETDAQ {
22 
24 
25 
26  /**
27  * Timeslice with random data.
28  */
30  public JTimesliceL0
31  {
32  /**
33  * Default constructor.
34  */
36  {}
37 
38 
39  /**
40  * Constructor.
41  *
42  * \param chronometer chronometer
43  * \param simbad detector simulator
44  */
45  JRandomTimeslice(const JDAQChronometer& chronometer,
46  const JDetectorSimulator& simbad)
47  {
48  using namespace JPP;
49 
50  setDAQChronometer(chronometer);
51 
52  if (simbad.hasK40Simulator() &&
53  simbad.hasPMTSimulator() &&
54  simbad.hasCLBSimulator()) {
55 
56  const double Tmin = getTimeSinceRTS(chronometer.getFrameIndex()); // [ns]
57 
58  const JTimeRange period(Tmin, Tmin + getFrameTime()); // [ns]
59 
60  JModuleData buffer;
61 
62  for (JDetector::const_iterator module = simbad->begin(); module != simbad->end(); ++module) {
63 
64  if (!module->empty()) {
65 
66  buffer.reset(module->size());
67 
68  simbad.generateHits(*module, period, buffer);
69 
70  this->push_back(JDAQSuperFrame(JDAQSuperFrameHeader(chronometer, module->getID())));
71 
72  simbad(*module, buffer, *(this->rbegin()));
73  }
74  }
75  }
76  }
77 
78 
79  /**
80  * Recycle time slice by randomly shuffling time intervals of data.
81  *
82  * Hits within one time interval are swapped with hits within another -randomly selected- time interval.
83  *
84  * Note that the time interval should be:
85  * - (much) larger than time differences of hits in L1 coincidence; and
86  * - (much) smaller than time covered by data (as per KM3NETDAQ::getFrameTime()).
87  *
88  * \param T_ns time interval
89  */
90  void recycle(const double T_ns)
91  {
92  using namespace std;
93  using namespace JPP;
94 
95  typedef JDAQHit::JTDC_t JTDC_t;
96 
97 
98  size_t N; // number of time intervals
99 
100  if (T_ns < 1.0)
101  N = (size_t) getFrameTime();
102  else if (T_ns < getFrameTime())
103  N = (size_t) (getFrameTime() / T_ns + 0.5);
104  else
105  N = 1;
106 
107  const JTDC_t T_step = (JTDC_t) (getFrameTime() / N); // TDC interval
108 
109 
110  typedef vector<JDAQHit> buffer_type;
111 
112  vector<size_t> index (N); // indices of time intervals for swapping
113  vector<buffer_type> buffer(N); // data per time interval
114 
115  for (iterator frame = this->begin(); frame != this->end(); ++frame) {
116 
117  if (!frame->empty()) {
118 
119 
120  for (size_t i = 0; i != N; ++i) {
121  index [i] = i;
122  buffer[i].clear();
123  }
124 
125 
126  vector<size_t> N_max(NUMBER_OF_PMTS, N); // number of time intervals per PMT
127  vector<JTDC_t> T_max(NUMBER_OF_PMTS, 0); // maximal time per PMT, e.g.\ due to high-rate veto
128 
129  for (JDAQSuperFrame::const_iterator hit = frame->begin(); hit != frame->end(); ++hit) {
130  T_max[hit->getPMT()] = hit->getT();
131  }
132 
133  for (int i = 0; i != NUMBER_OF_PMTS; ++i) {
134 
135  const size_t N_i = (T_max[i] + T_step) / T_step;
136 
137  if (N_i < N_max[i]) {
138  N_max[i] = N_i;
139  }
140  }
141 
142  sort(N_max.begin(), N_max.end()); // swaps by allowed time range
143 
144  for (size_t i = 0, L = 0; i != N_max.size(); ++i) {
145 
146  size_t M = N_max[i];
147 
148  if (L != M) {
149  for (size_t i = M - L - 1; i != 0; --i) {
150  std::swap(index[L + i], index[L + gRandom->Integer(i+1)]);
151  }
152  }
153 
154  L = M;
155  }
156 
157  const size_t M = N_max[NUMBER_OF_PMTS - 1];
158 
159 
160  // Wall street shuflle
161 
162  JDAQSuperFrame::const_iterator hit = frame->begin();
163 
164  for (size_t in = 0; in != M; ++in) {
165 
166  const size_t out = index[in];
167  buffer_type& zbuf = buffer[out];
168 
169  const JTDC_t T_in = in * T_step;
170  const JTDC_t T_out = out * T_step;
171 
172  for ( ; hit != frame->end() && hit->getT() < T_in + T_step; ++hit) {
173  zbuf.push_back(JDAQHit(hit->getPMT(), (hit->getT() - T_in) + T_out, hit->getToT()));
174  }
175  }
176 
177 
178  // copy back data
179 
180  for (size_t i = 0, number_of_hits = 0; i != M; ++i) {
181 
182  memcpy(frame->data() + number_of_hits, buffer[i].data(), buffer[i].size() * sizeof(JDAQHit));
183 
184  number_of_hits += buffer[i].size();
185  }
186  }
187  }
188  }
189  };
190 }
191 
192 #endif
193 
void reset(size_t size)
Reset buffers.
Data structure for PMT data corresponding to a detector module.
JRandomTimeslice()
Default constructor.
Auxiliary class for TDC constraints.
Definition: JTDC_t.hh:37
then JShowerPostfit f $INPUT_FILE o $OUTPUT_FILE N
JTDC_t getT() const
Get time.
Definition: JDAQHit.hh:86
JPMT_t getPMT() const
Get PMT.
Definition: JDAQHit.hh:75
unsigned int JTDC_t
leading edge [ns]
Definition: JDAQHit.hh:39
int getFrameIndex() const
Get frame index.
virtual void generateHits(const JModule &module, const JTimeRange &period, JModuleData &output) const
Generate hits.
JTOT_t getToT() const
Get time-over-threshold.
Definition: JDAQHit.hh:97
bool hasPMTSimulator() const
Check availability of PMT simulator.
Hit data structure.
Definition: JDAQHit.hh:34
double getFrameTime()
Get frame time duration.
Definition: JDAQClock.hh:162
void setDAQChronometer(const JDAQChronometer &chronometer)
Set DAQ chronometer.
bool hasK40Simulator() const
Check availability of K40 simulator.
bool hasCLBSimulator() const
Check availability of CLB simulator.
double getTimeSinceRTS(const int frame_index)
Get time in ns since last RTS for a given frame index.
Definition: JDAQClock.hh:263
JRandomTimeslice(const JDAQChronometer &chronometer, const JDetectorSimulator &simbad)
Constructor.
Base class class for generation of time slice data.
Definition: JTimesliceL0.hh:18
Auxiliaries for creation of time slice data.
KM3NeT DAQ constants, bit handling, etc.
static const int NUMBER_OF_PMTS
Total number of PMTs in module.
Definition: JDAQ.hh:26
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
Timeslice with random data.
Data frame of one optical module.
void recycle(const double T_ns)
Recycle time slice by randomly shuffling time intervals of data.