Project

General

Profile

dump_gsimple_file.C

script to dump GSimpleNtpFlux file entries - Robert Hatcher, 07/25/2012 06:15 PM

 
1
#include "Numerical/RandomGen.h"
2
#include "FluxDrivers/GSimpleNtpFlux.h"
3
#include "Utils/UnitUtils.h"
4

    
5
#include "TSystem.h"
6
#include "TStopwatch.h"
7
#include "TLorentzVector.h"
8
#include "TNtuple.h"
9
#include "TFile.h"
10
#include "TChain.h"
11

    
12
#include <iostream>
13
#include <iomanip>
14
#include <string>
15
#include <sstream>
16
#include <set>
17

    
18
using namespace std;
19
using namespace genie;
20
using namespace genie::flux;
21

    
22
//========================================================================
23
// main routine
24

    
25
// assumes a single input file (meta-data handling)
26

    
27

    
28
void dump_gsimple_file(string fname="gsimple_output.root",
29
                       long int nentries=10)
30
{
31
    string fluxfname(gSystem->ExpandPathName(fname.c_str()));
32
  flux::GSimpleNtpFlux* gsimplein = new GSimpleNtpFlux();
33
  gsimplein->LoadBeamSimData(fluxfname,"<<no-offset-index>>");
34
  gsimplein->SetEntryReuse(1); // don't reuse entries when dumping
35
  gsimplein->GenerateWeighted(true); // don't deweight, we want to see file entries
36
  gsimplein->SetUpstreamZ(-3e38);  // leave ray on flux window
37

    
38
  const string sepline = 
39
    "========================================================================";
40

    
41
  cout << sepline << endl << flush;
42
  gsimplein->PrintConfig();
43
  cout << sepline << endl << flush;
44

    
45
  GFluxI* fdriver = dynamic_cast<GFluxI*>(gsimplein);
46

    
47
  // if unspecified do all the entries in the input file
48
  if (nentries <= 0) {
49
    // nentries = 2147483647;
50
    nentries = gsimplein->GetFluxTChain()->GetEntries();
51
  }
52

    
53
  UInt_t last_metakey = 0;
54
  for ( long int ientry = 0; ientry < nentries; ++ientry ) {
55
    fdriver->GenerateNext();
56

    
57
    cout << *(gsimplein->GetCurrentEntry())
58
         << *(gsimplein->GetCurrentAux())
59
         << *(gsimplein->GetCurrentNuMI())
60
         << endl << flush;
61

    
62
    const genie::flux::GSimpleNtpMeta* fmeta_in = gsimplein->GetCurrentMeta();
63
    if ( last_metakey != fmeta_in->metakey ) {
64
      cout << *fmeta_in << endl << flush;
65
      last_metakey = fmeta_in->metakey;
66
    }
67
 
68
  }
69

    
70
  cout << "=========================== Complete " << endl;
71

    
72
  cout << "Generated/Dumped " << nentries << " entries"
73
       << endl
74
       << gsimplein->UsedPOTs() << " POTs " 
75
       << ", pulled NFluxNeutrinos " << gsimplein->NFluxNeutrinos()
76
       << endl;
77
  
78
}
79