Jpp test-rotations-old
the software that should make you happy
Loading...
Searching...
No Matches
JSUPERNOVA::JLightCurveBackgroundGenerator Class Reference

Class to emulate L0 background for an arbitrarily sized detector. More...

#include <JLightCurveBackgroundGenerator.hh>

Public Member Functions

 JLightCurveBackgroundGenerator (TFile *in, bool twoDim=false)
 Default constructor.
 
 ~JLightCurveBackgroundGenerator ()
 Destructor.
 
void configureTimeWindow (const int T_ms)
 Configure the duration of an output sample.
 
void setSeed (const UInt_t uSeed=0)
 Set TRandom seed.
 
void configureRatio (const int inputNumberOfLines, const int outputNumberOfLines)
 Configure generation ratio.
 
bg_type generate ()
 Generate sample of L0 background L0 data are randomly sampled from a single L0 dataset.
 
bg_type generate_poisson (const int domRate_Hz=500)
 Generate pure poissionian background.
 
bg_type generate_shuffled (bool randomizeRun=false)
 Generate sample of L0 background The sampling of the L0 data is not sequential but random within the L0 dataset.
 
h2d_tbuild_H2D (const int k, const int start)
 Generate 2D histogram.
 
vector< double > build_NC (const int k, const int start)
 Build NC sequence.
 
h2d_tgenerate_H2D (int rb=1)
 Generate 2D sample.
 
bg_type generate_fitted (int rb=1)
 Generate fitted L1 sample.
 

Static Public Member Functions

static vector< double > rebin (const vector< double > &in, const int rb)
 Rebin vector.
 
static bg_type fit_H2D (h2d_t *in)
 Fit 2D sample.
 

Private Attributes

vector< vector< double > > rt
 
vector< vector< double > > nc
 
vector< h2d_t * > rt2d
 
int TWindow_ms
 
int genRatio
 
int outputSize
 
double genScale
 
TRandom * rnd
 

Detailed Description

Class to emulate L0 background for an arbitrarily sized detector.

Definition at line 59 of file JLightCurveBackgroundGenerator.hh.

Constructor & Destructor Documentation

◆ JLightCurveBackgroundGenerator()

JSUPERNOVA::JLightCurveBackgroundGenerator::JLightCurveBackgroundGenerator ( TFile * in,
bool twoDim = false )
inline

Default constructor.

The input file contains a sequence of histograms with sampled L0 data; Each histogram is loaded as a vector, which will be considered as an individual L0 dataset.

Parameters
ininput file from JRipple
twoDimmake 2D histogram with distribution of hit time differences vs time

Definition at line 85 of file JLightCurveBackgroundGenerator.hh.

85 {
86
87 TWindow_ms = 100;
88 genRatio = 115;
89 genScale = 1;
90
91 TString rt_tag = ".RT_DET_SUM";
92 TString nc_tag = ".NC_DET_SUM";
93 TString rt2d_tag = ".RT2D_DET";
94
95 rnd = new TRandom();
96
97 TIter iter(in->GetListOfKeys());
98
99 for (TKey* key; (key = (TKey*) iter.Next()) != NULL; ) {
100
101 const TString tag(key->GetName());
102
103 // match
104 if (tag.EndsWith(rt_tag)) {
105
106 TString rt_hist_tag = tag;
107 TString run_tag = TString(tag, tag.Length() - rt_tag.Length());
108
109 TString nc_hist_tag = run_tag + nc_tag;
110
111 TH1D* RT = (TH1D*) in->Get(rt_hist_tag);
112 TH1D* NC = (TH1D*) in->Get(nc_hist_tag);
113
114 vector<double> rt_buf = loadHistogram(RT);
115 vector<double> nc_buf = loadHistogram(NC);
116
117 int n = nc_buf.size();
118
119 rt.push_back( vector<double>() );
120 nc.push_back( vector<double>() );
121
122 // build parallel vectors
123 // 1 NC bin <=> 100 RT bins
124 for (int i = 0; i < n; i++) {
125 // empty slices are discarded
126 if (nc_buf[i] > 0) {
127 nc.back().push_back( nc_buf[i] );
128 for (int j = 100 * i; j < 100 * (i + 1); j++) {
129 rt.back().push_back( rt_buf[j] );
130 }
131 }
132 }
133
134 if (twoDim) {
135 TString rt2d_hist_tag = run_tag + rt2d_tag;
136 h2d_t* RT2D = (h2d_t*) in->Get(rt2d_hist_tag);
137
138 if (RT2D != NULL) {
139 RT2D->SetDirectory(0);
140 rt2d.push_back(RT2D);
141 }
142 }
143
144 }
145
146 }
147
148 }
vector< double > loadHistogram(TH1D *in)
Load histogram values to vector, each bin is converted to an element.
TH2F h2d_t
Definition JRipple.hh:7
const int n
Definition JPolint.hh:791
int j
Definition JPolint.hh:801

◆ ~JLightCurveBackgroundGenerator()

JSUPERNOVA::JLightCurveBackgroundGenerator::~JLightCurveBackgroundGenerator ( )
inline

Destructor.

Definition at line 153 of file JLightCurveBackgroundGenerator.hh.

153 {
154 delete rnd;
155 }

Member Function Documentation

◆ configureTimeWindow()

void JSUPERNOVA::JLightCurveBackgroundGenerator::configureTimeWindow ( const int T_ms)
inline

Configure the duration of an output sample.

Parameters
T_mssample duration in ms

Definition at line 162 of file JLightCurveBackgroundGenerator.hh.

162 {
163 TWindow_ms = T_ms;
164 }

◆ setSeed()

void JSUPERNOVA::JLightCurveBackgroundGenerator::setSeed ( const UInt_t uSeed = 0)
inline

Set TRandom seed.

Parameters
uSeedseed

Definition at line 171 of file JLightCurveBackgroundGenerator.hh.

171 {
172 rnd->SetSeed(uSeed);
173 }

◆ configureRatio()

void JSUPERNOVA::JLightCurveBackgroundGenerator::configureRatio ( const int inputNumberOfLines,
const int outputNumberOfLines )
inline

Configure generation ratio.

Parameters
inputNumberOfLinesnumber of lines of input detector
outputNumberOfLinesnumber of lines of output detector

Definition at line 181 of file JLightCurveBackgroundGenerator.hh.

181 {
182 outputSize = inputNumberOfLines;
183 // ratio is integer by design
184 genRatio = outputNumberOfLines / inputNumberOfLines;
185 // scaling factor to compensate the ratio remainder
186 genScale = outputNumberOfLines / (1.0 * inputNumberOfLines * genRatio);
187 }

◆ generate()

bg_type JSUPERNOVA::JLightCurveBackgroundGenerator::generate ( )
inline

Generate sample of L0 background L0 data are randomly sampled from a single L0 dataset.

Returns
2D vector with hit count and number of active channels as a function of time

Definition at line 196 of file JLightCurveBackgroundGenerator.hh.

196 {
197
198 bg_type out(2, vector<double>(TWindow_ms));
199
200 int k = rnd->Integer(rt.size());
201
202 int limit = rt[k].size() - (genRatio * TWindow_ms);
203
204 int start = rnd->Integer(limit);
205
206 for (int i = 0; i < genRatio; i++) {
207
208 for (int j = 0; j < TWindow_ms; j++) {
209
210 int bin = start + (i * TWindow_ms) + j;
211
212 out[0][j] += genScale * rt[k].at(bin );
213 out[1][j] += genScale * nc[k].at(bin / 100);
214
215 }
216
217 }
218
219 return out;
220 }
vector< vector< double > > bg_type

◆ generate_poisson()

bg_type JSUPERNOVA::JLightCurveBackgroundGenerator::generate_poisson ( const int domRate_Hz = 500)
inline

Generate pure poissionian background.

Definition at line 225 of file JLightCurveBackgroundGenerator.hh.

225 {
226
227 bg_type out(2, vector<double>(TWindow_ms));
228
229 double mu = outputSize * 18 * domRate_Hz * 1.0e-3;
230
231 for (unsigned i = 0; i < out[0].size(); i++) {
232
233 double count = rnd->Poisson(mu);
234
235 out[0][i] = count;
236 out[1][i] = outputSize * 18 * 31;
237
238 }
239
240 return out;
241
242 }

◆ generate_shuffled()

bg_type JSUPERNOVA::JLightCurveBackgroundGenerator::generate_shuffled ( bool randomizeRun = false)
inline

Generate sample of L0 background The sampling of the L0 data is not sequential but random within the L0 dataset.

Returns
2D vector with hit count and number of active channels as a function of time

Definition at line 250 of file JLightCurveBackgroundGenerator.hh.

250 {
251
252 bg_type out(2, vector<double>(TWindow_ms));
253
254 int nRuns = rt.size();
255
256 int k = rnd->Integer(nRuns);
257
258 // a series n = 'genRatio' elements is sampled randomly within the dataset
259 // all the elements start with the same offset to the timeslice beginning
260
261 int offset = rnd->Integer(100);
262
263 for (int i = 0; i < genRatio; i++) {
264
265 int start = 100 * rnd->Integer(nc[k].size() - ceil(TWindow_ms / 100.0)) + offset;
266
267 for (int j = 0; j < TWindow_ms; j++) {
268
269 int bin = start + j;
270
271 out[0][j] += genScale * rt[k].at(bin );
272 out[1][j] += genScale * nc[k].at(bin / 100);
273
274 }
275
276 if (randomizeRun) {
277 k = rnd->Integer(nRuns);
278 }
279
280 }
281
282 return out;
283 }

◆ build_H2D()

h2d_t * JSUPERNOVA::JLightCurveBackgroundGenerator::build_H2D ( const int k,
const int start )
inline

Generate 2D histogram.

Definition at line 289 of file JLightCurveBackgroundGenerator.hh.

289 {
290
291 h2d_t* proto = rt2d[k];
292
293 int ny = proto->GetNbinsY();
294 int ymin = proto->GetYaxis()->GetXmin();
295 int ymax = proto->GetYaxis()->GetXmax();
296
297 h2d_t* out = new h2d_t("sample", NULL, TWindow_ms, 0, TWindow_ms, ny, ymin, ymax);
298
299 for (int i = 0; i < genRatio; i++) {
300
301 for (int j = 0; j < TWindow_ms; j++) {
302
303 int xbin = 1 + start + (i * TWindow_ms) + j;
304
305 for (int ybin = 1; ybin <= rt2d[k]->GetNbinsY(); ybin++) {
306
307 double val = rt2d[k]->GetBinContent(xbin, ybin);
308
309 double bcx = ((TAxis*)out->GetXaxis())->GetBinCenter(j + 1);
310 double bcy = ((TAxis*)out->GetYaxis())->GetBinCenter(ybin);
311
312 out->Fill(bcx, bcy, val);
313
314 }
315
316 }
317
318 }
319
320 return out;
321
322 }

◆ build_NC()

vector< double > JSUPERNOVA::JLightCurveBackgroundGenerator::build_NC ( const int k,
const int start )
inline

Build NC sequence.

Definition at line 328 of file JLightCurveBackgroundGenerator.hh.

328 {
329
330 vector<double> out(TWindow_ms);
331
332 for (int i = 0; i < genRatio; i++) {
333
334 for (int j = 0; j < TWindow_ms; j++) {
335
336 int bin = start + (i * TWindow_ms) + j;
337
338 out[j] += genScale * nc[k].at(bin / 100);
339
340 }
341
342 }
343
344 return out;
345 }

◆ rebin()

static vector< double > JSUPERNOVA::JLightCurveBackgroundGenerator::rebin ( const vector< double > & in,
const int rb )
inlinestatic

Rebin vector.

Definition at line 350 of file JLightCurveBackgroundGenerator.hh.

350 {
351
352 const int n = in.size();
353
354 vector<double> out(n/rb);
355
356 for (int i = 0; i < n; i++) {
357 out[i/rb] += in[i];
358 }
359
360 return out;
361
362 }

◆ generate_H2D()

h2d_t * JSUPERNOVA::JLightCurveBackgroundGenerator::generate_H2D ( int rb = 1)
inline

Generate 2D sample.

Definition at line 367 of file JLightCurveBackgroundGenerator.hh.

367 {
368
369 h2d_t* out = NULL;
370
371 if (rt2d.size() > 0) {
372
373 int k = rnd->Integer(rt2d.size());
374
375 int limit = rt2d[k]->GetNbinsX() - (genRatio * TWindow_ms);
376
377 int start = rnd->Integer(limit);
378
379 out = build_H2D(k, start);
380
381 }
382
383 if (rb != 1) { out->RebinX(rb); }
384
385 return out;
386
387 }
h2d_t * build_H2D(const int k, const int start)
Generate 2D histogram.

◆ fit_H2D()

static bg_type JSUPERNOVA::JLightCurveBackgroundGenerator::fit_H2D ( h2d_t * in)
inlinestatic

Fit 2D sample.

Definition at line 393 of file JLightCurveBackgroundGenerator.hh.

393 {
394
395 int n = in->GetNbinsX();
396
397 // 0 -> signal
398 // 1 -> fitted background
399 // 2 -> background fit error
400
401 bg_type out(3, vector<double>(n));
402
403 for (int bin = 1; bin <= n; bin++) {
404
405 double sg = 0, bg = 0, er = 0;
406
407 TH1D* h = in->ProjectionY("timeBin", bin, bin);
408
409 if (h->GetEntries() > 0) {
410
411 TF1 f("f", "[0]*exp(-0.5*(x-[1])*(x-[1])/([2]*[2]))/(TMath::Sqrt(2*TMath::Pi())*[2]) + [3]");
412
413 f.SetParameter(0, h->GetMaximum());
414 f.SetParameter(1, h->GetMean());
415 f.SetParameter(2, h->GetRMS() * 0.25);
416 f.SetParameter(3, h->GetMinimum());
417
418 h->Fit(&f, "Q", "same");
419
420 bg = f.GetParameter(3) * h->GetNbinsX();
421 sg = h->GetSumOfWeights() - bg;
422 er = f.GetParError(3) * h->GetNbinsX();
423
424 }
425
426 int i = bin - 1;
427
428 out[0][i] = sg;
429 out[1][i] = bg;
430 out[2][i] = er;
431
432 delete h;
433
434 }
435
436 return out;
437
438 }

◆ generate_fitted()

bg_type JSUPERNOVA::JLightCurveBackgroundGenerator::generate_fitted ( int rb = 1)
inline

Generate fitted L1 sample.

Definition at line 445 of file JLightCurveBackgroundGenerator.hh.

445 {
446
447 bg_type out;
448
449 int k = rnd->Integer(rt.size());
450
451 int limit = rt[k].size() - (genRatio * TWindow_ms);
452
453 int start = rnd->Integer(limit);
454
455 h2d_t* sample = build_H2D(k, start);
456
457 if (rb != 1) { sample->RebinX(rb); }
458
459 bg_type fit = fit_H2D(sample);
460
461 out.push_back(fit[0]);
462
463 delete sample;
464
465 if (rb == 1) {
466 out.push_back(build_NC(k, start));
467 } else {
468 out.push_back(rebin(build_NC(k, start), rb));
469 }
470
471 out.push_back(fit[1]);
472 out.push_back(fit[2]);
473
474 return out;
475
476 }
vector< double > build_NC(const int k, const int start)
Build NC sequence.
static vector< double > rebin(const vector< double > &in, const int rb)
Rebin vector.

Member Data Documentation

◆ rt

vector<vector<double> > JSUPERNOVA::JLightCurveBackgroundGenerator::rt
private

Definition at line 61 of file JLightCurveBackgroundGenerator.hh.

◆ nc

vector<vector<double> > JSUPERNOVA::JLightCurveBackgroundGenerator::nc
private

Definition at line 62 of file JLightCurveBackgroundGenerator.hh.

◆ rt2d

vector<h2d_t*> JSUPERNOVA::JLightCurveBackgroundGenerator::rt2d
private

Definition at line 64 of file JLightCurveBackgroundGenerator.hh.

◆ TWindow_ms

int JSUPERNOVA::JLightCurveBackgroundGenerator::TWindow_ms
private

Definition at line 66 of file JLightCurveBackgroundGenerator.hh.

◆ genRatio

int JSUPERNOVA::JLightCurveBackgroundGenerator::genRatio
private

Definition at line 67 of file JLightCurveBackgroundGenerator.hh.

◆ outputSize

int JSUPERNOVA::JLightCurveBackgroundGenerator::outputSize
private

Definition at line 68 of file JLightCurveBackgroundGenerator.hh.

◆ genScale

double JSUPERNOVA::JLightCurveBackgroundGenerator::genScale
private

Definition at line 70 of file JLightCurveBackgroundGenerator.hh.

◆ rnd

TRandom* JSUPERNOVA::JLightCurveBackgroundGenerator::rnd
private

Definition at line 72 of file JLightCurveBackgroundGenerator.hh.


The documentation for this class was generated from the following file: