star-travex
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
StiScanRootFile.cxx
Go to the documentation of this file.
1 #include <cmath>
2 #include <boost/filesystem.hpp>
3 
4 #include "StiScanRootFile.h"
5 
6 #include "TRandom.h"
7 
12 
13 
14 StiScanRootFile::StiScanRootFile(StiScanPrgOptions& prgOpts, Option_t *option, const char *ftitle, Int_t compress) :
15  tvx::RootFile(prgOpts, option, ftitle, compress),
16  fPrgOptions(prgOpts)
17 {
18  // Find ranges (\todo if requested by the user)
20  Info("StiScanRootFile", "Find auto range. Loop over tree/chain...");
21  FindAutoRange();
22  }
23 
24  Add( new StiScanHistsByVolume(fPrgOptions, "sti_vol", this) );
25  Add( new StiScanHistContainer(fPrgOptions, "sti_trk", this) );
26  Add( new StiScanHistContainer(fPrgOptions, "gea", this) ); // Will create integral projections from 2D histograms
27 }
28 
29 
31 {
32  const std::set<std::string> *volumeList = &fPrgOptions.GetVolumeList();
33 
34  TChain *stiChain = fPrgOptions.GetStiTChain();
35 
36  StiScanEvent *stiEvent = new StiScanEvent();
37  stiChain->SetBranchAddress("e.", &stiEvent);
38  stiChain->SetBranchStatus("e.TStiEvent.*", false);
39  stiChain->SetBranchStatus("e.TStiEvent.fTStiKalmanTracks*", true);
40 
41  double nodeRMax = 0;
42  double nodeZMin = DBL_MAX;
43  double nodeZMax = -DBL_MAX;
44 
45  TRandom myRandom;
46 
47  int nTreeEvents = stiChain->GetEntries();
48 
49  Info("FindAutoRange", "Found tree/chain with N entries: %d", nTreeEvents);
50 
51  for (int iEvent = 1; iEvent <= nTreeEvents; iEvent++)
52  {
53  if ( nTreeEvents >= 10 && iEvent % int(nTreeEvents*0.1) == 0 )
54  Info("FindAutoRange", "Analyzing event %d", iEvent);
55 
56  if (myRandom.Rndm() > fPrgOptions.GetSparsity()) continue;
57 
58  stiChain->GetEntry(iEvent-1);
59 
60  for (const auto& kalmTrack : stiEvent->GetTStiKalmanTracks())
61  {
62  for (const auto& kalmNode : kalmTrack.GetNodes())
63  {
64  if (volumeList && volumeList->size() && !kalmNode.MatchedVolName(*volumeList) )
65  continue;
66 
67  if (kalmNode.GetNodeMaterialDensity() <= 0)
68  continue;
69 
70  double node_z = kalmNode.GetPosition().Z();
71 
72  if (node_z < nodeZMin) nodeZMin = node_z;
73  if (node_z > nodeZMax) nodeZMax = node_z;
74 
75  // Find maximum radius
76  double node_r = (float) kalmNode.GetPosition().Perp();
77 
78  if (node_r > nodeRMax) nodeRMax = node_r;
79  }
80  }
81  }
82 
83  // In case no entries were found in the tree set some reasonable values
84  if (nodeRMax == 0 && nodeZMin == DBL_MAX && nodeZMax == -DBL_MAX) {
85  nodeRMax = 1;
86  nodeZMin = 0;
87  nodeZMax = 1;
88  }
89 
90  Info("FindAutoRange", "Updated nodeZMin, nodeZMax, nodeRMax: %g, %g, %g", nodeZMin, nodeZMax, nodeRMax);
91 
92  fPrgOptions.SetHistZRange(nodeZMin, nodeZMax);
93  fPrgOptions.SetHistRMax(nodeRMax);
94 
95  delete stiEvent;
96 }
97 
98 
99 void StiScanRootFile::FillHists(const StiScanEvent &stiEvent, const std::set<std::string> *volumeList)
100 {
101  hc<StiScanHistContainer>("sti_vol")->FillHists(stiEvent, volumeList);
102  hc<StiScanHistContainer>("sti_trk")->FillHists(stiEvent, volumeList);
103 }
104 
105 
106 void StiScanRootFile::FillHists(const TGeaEvent &geaEvent, const std::set<std::string> *volumeList)
107 {
108  hc<StiScanHistContainer>("gea")->FillHists(geaEvent, volumeList);
109 }
110 
111 
113 {
114  tvx::RootFile::Finalize();
115 
116  StiScanRatiosHistContainer *ratios = new StiScanRatiosHistContainer("sti_gea_ratio", this);
117  Add(ratios);
118 
119  const tvx::HistContainer& gea = *hc("gea");
120  const tvx::HistContainer& sti_trk = *hc("sti_trk");
121 
122  const TH1* gea_eloss = &gea["hELossVsPhiVsRVsZ_pyx"];
123  const TH1* sti_eloss = &sti_trk["hELossVsPhiVsRVsZ_pyx"];
124  ratios->CreateRatioHist(sti_eloss, gea_eloss);
125 
126  gea_eloss = &gea["hELossVsPhiVsRVsZ_pyx_px"];
127  sti_eloss = &sti_trk["hELossVsPhiVsRVsZ_pyx_px"];
128  ratios->CreateRatioHist(sti_eloss, gea_eloss);
129 
130  gea_eloss = &gea["hELossVsXVsYVsZ_pyx"];
131  sti_eloss = &sti_trk["hELossVsXVsYVsZ_pyx"];
132  ratios->CreateRatioHist(sti_eloss, gea_eloss);
133 }
bool DoAutoHistRange() const
void FindAutoRange() const
StiScanPrgOptions & fPrgOptions
Shadow base class reference.
void SetHistRMax(double maxR)
void CreateRatioHist(const TH1 *hNumer, const TH1 *hDenom)
void FillHists(const StiScanEvent &stiEvent, const std::set< std::string > *volumeList=0)
void SetHistZRange(double minZ, double maxZ)
Processes and controls user options provided in the command line.
const std::vector< TStiKalmanTrack > & GetTStiKalmanTracks() const
Definition: TStiEvent.h:38
StiScanRootFile(StiScanPrgOptions &prgOpts, Option_t *option="", const char *ftitle="", Int_t compress=1)
const std::set< std::string > & GetVolumeList() const