Project

General

Profile

Example Analysis Code » cpp_analysis_example.cc

Example C++ analysis - Dobbs, Adam, 08 May 2015 17:32

 
1
#include <iostream>
2
#include <string>
3
#include <vector>
4

    
5
#include "TCanvas.h"
6
#include "TH1D.h"
7
#include "TApplication.h"
8

    
9
#include "src/common_cpp/JsonCppStreamer/IRStream.hh"
10
#include "src/common_cpp/DataStructure/Spill.hh"
11
#include "src/common_cpp/DataStructure/Data.hh"
12
#include "src/common_cpp/DataStructure/ReconEvent.hh"
13
#include "src/common_cpp/DataStructure/SciFiEvent.hh"
14
#include "src/common_cpp/DataStructure/SciFiSpacePoint.hh"
15

    
16
/** Access Tracker data using ROOT */
17

    
18
int main(int argc, char *argv[]) {
19
  // First argument to code should be the input ROOT file name
20
  std::string filename = std::string(argv[1]);
21

    
22
  // Set up ROOT app, input file, and MAUS data class
23
  TApplication theApp("App", &argc, argv);
24
  std::cout << "Opening file " << filename << std::endl;
25
  irstream infile(filename.c_str(), "Spill");
26
  MAUS::Data data;
27

    
28
  // Create the histogram
29
  TH1D* h1 = new TH1D("npe", "SciFi Spacepoint NPE", 100, 0, 150);
30
  h1->GetXaxis()->SetTitle("Number of photoelectrons");
31
  h1->SetLineColor(kBlue);
32

    
33
  // Loop over all spills
34
  while ( infile >> readEvent != NULL ) {
35
    infile >> branchName("data") >> data;
36
    MAUS::Spill* spill = data.GetSpill();
37
    if (spill == NULL || spill->GetDaqEventType() != "physics_event")
38
      continue;
39
    std::cout << "Spill: " << spill->GetSpillNumber() << "\n";
40
    std::vector<MAUS::ReconEvent*>* revts = spill->GetReconEvents();
41

    
42
    // Loop over recon events in spill
43
    for ( size_t i = 0; i < revts->size(); ++i ) {
44
      if ( !revts->at(i) ) continue;
45
      MAUS::SciFiEvent* sfevt = revts->at(i)->GetSciFiEvent();
46

    
47
      // Loop over spacepoints
48
      std::vector<MAUS::SciFiSpacePoint*> spnts = sfevt->spacepoints();
49
      std::vector<MAUS::SciFiSpacePoint*>::iterator spnt;
50
      for ( spnt = spnts.begin(); spnt != spnts.end(); ++spnt ) {
51
        h1->Fill((*spnt)->get_npe());
52
      }
53
    } // ~Loop over Recon events
54
  } // ~Loop over all spills
55

    
56
  // Draw the histogram
57
  TCanvas * c1 = new TCanvas("c1", "SciFi Spacepoint NPE", 900, 600);
58
  c1->cd();
59
  h1->Draw();
60
  theApp.Run();
61
}
62

    
(2-2/5)