BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
rebdsimCombine.cc
Go to the documentation of this file.
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2022.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
22#include "FileMapper.hh"
23#include "Header.hh"
24#include "HistogramAccumulatorMerge.hh"
25#include "HistogramAccumulatorSum.hh"
26#include "RBDSException.hh"
27
28#include "BDSOutputROOTEventHeader.hh"
29
30#include "TDirectory.h"
31#include "TFile.h"
32#include "TH1.h"
33#include "TH2.h"
34#include "TH3.h"
35#include "BDSBH4D.hh"
36#include "TTree.h"
37
38#include <exception>
39#include <iostream>
40#include <set>
41#include <string>
42#include <vector>
43
44int main(int argc, char* argv[])
45{
46 if (argc < 3)
47 {
48 std::cout << "usage: rebdsimCombine result.root file1.root file2.root ..." << std::endl;
49 return 1;
50 }
51
52 // build input file list
53 std::vector<std::string> inputFiles;
54 for (int i = 2; i < argc; ++i)
55 {inputFiles.emplace_back(std::string(argv[i]));}
56
57 // checks
58 if (inputFiles.size() == 1)
59 {
60 if (inputFiles[0].find('*') != std::string::npos) // glob didn't expand in shell - infer this
61 {std::cout << "Glob with * did not match any files" << std::endl;}
62 else
63 {std::cout << "Only one input file provided \"" << inputFiles[0] << "\" - no point." << std::endl;}
64 return 1;
65 }
66
67 std::set<std::string> inputFilesSet = {inputFiles.begin(), inputFiles.end()};
68 if (inputFiles.size() > inputFilesSet.size())
69 {std::cout << "Warning: at least 1 duplicate name in list of files provided to combine." << std::endl;}
70
71 std::string outputFile = std::string(argv[1]);
72 if (outputFile.find('*') != std::string::npos)
73 {
74 std::cerr << "First argument for output file \"" << outputFile << "\" contains an *." << std::endl;
75 std::cerr << "Should only be a singular file - check order of arguments." << std::endl;
76 return 1;
77 }
78
79 // output file must be opened before histograms are created because root does
80 // everything statically behind the scenes
81 TFile* output = new TFile(outputFile.c_str(), "RECREATE");
82
83 // add header for file type and version details
84 output->cd();
86 headerOut->Fill(std::vector<std::string>(), inputFiles); // updates time stamp
87 headerOut->SetFileType("REBDSIMCOMBINE");
88 TTree* headerTree = new TTree("Header", "REBDSIM Header");
89 headerTree->Branch("Header.", "BDSOutputROOTEventHeader", headerOut);
90
91 // ensure new histograms are written to file
92 TH1::AddDirectory(true);
93 TH2::AddDirectory(true);
94 TH3::AddDirectory(true);
95
96 TFile* f = nullptr; // temporary variable
97
98 // initialise file map
99 try
100 {f = new TFile(inputFiles[0].c_str(), "READ");}
101 catch (const std::exception& e)
102 {std::cerr << e.what() << std::endl; return 1;}
103 HistogramMap* histMap = nullptr;
104 try
105 {histMap = new HistogramMap(f, output);} // map out first file
106 catch (const RBDSException& error)
107 {std::cerr << error.what() << std::endl; return 1;}
108 catch (const std::exception& error)
109 {std::cerr << error.what() << std::endl; return 1;}
110
111 // copy the model tree over if it exists - expect the name to be ModelTree
112 TTree* oldModelTree = dynamic_cast<TTree*>(f->Get("ModelTree"));
113 if (!oldModelTree)
114 {oldModelTree = dynamic_cast<TTree*>(f->Get("Model"));}
115 if (oldModelTree)
116 {// TChain can be valid but TTree might not be in corrupt / bad file
117 output->cd();
118 auto newTree = oldModelTree->CloneTree();
119 newTree->SetName("ModelTree");
120 newTree->Write("", TObject::kOverwrite);
121 }
122
123 f->Close();
124 delete f;
125
126 std::vector<RBDS::HistogramPath> histograms = histMap->Histograms();
127
128 unsigned long long int nOriginalEvents = 0;
129
130 std::cout << "Combination of " << inputFiles.size() << " files beginning" << std::endl;
131 // loop over files and accumulate
132 for (const auto& file : inputFiles)
133 {
134 f = new TFile(file.c_str());
136 {
137 std::cout << "Accumulating> " << file << std::endl;
138 for (const auto& hist : histograms)
139 {
140 std::string histPath = hist.path + hist.name; // histPath has trailing '/'
141
142 TH1* h = dynamic_cast<TH1*>(f->Get(histPath.c_str()));
143
144 if (!h)
145 {RBDS::WarningMissingHistogram(histPath, file); continue;}
146 hist.accumulator->Accumulate(h);
147 }
148
149 Header* h = new Header();
150 TTree* ht = (TTree*)f->Get("Header");
151 h->SetBranchAddress(ht);
152 ht->GetEntry(0);
153 nOriginalEvents += h->header->nOriginalEvents;
154 delete h;
155 }
156 else
157 {std::cout << "Skipping " << file << " as not a rebdsim output file" << std::endl;}
158 f->Close();
159 delete f;
160 }
161
162 // terminate and write output
163 for (const auto& hist : histograms)
164 {
165 TH1* result = hist.accumulator->Terminate();
166 result->SetDirectory(hist.outputDir);
167 hist.outputDir->Add(result);
168 delete hist.accumulator; // this removes temporary histograms from the file
169 }
170
171 headerOut->nOriginalEvents = nOriginalEvents;
172 headerTree->Fill();
173
174 output->Write(nullptr,TObject::kOverwrite);
175 output->Close();
176 delete output;
177 std::cout << "Combined result of " << inputFiles.size() << " files written to: " << outputFile << std::endl;
178
179 return 0;
180}
Information about the software and the file.
void SetFileType(std::string fileTypeIn)
Update the file type.
unsigned long long int nOriginalEvents
Number of original events if skimmed.
void Fill(const std::vector< std::string > &analysedFilesIn=std::vector< std::string >(), const std::vector< std::string > &combinedFilesIn=std::vector< std::string >())
Header loader.
Definition: Header.hh:34
BDSOutputROOTEventHeader * header
Member that ROOT can map file data to locally.
Definition: Header.hh:44
void SetBranchAddress(TTree *t)
Set the branch addresses to address the contents of the file.
Definition: Header.cc:41
Class to map a rebdsim file structure and create duplicate in output.
Definition: FileMapper.hh:95
const std::vector< RBDS::HistogramPath > & Histograms() const
Access full vector of histograms.
Definition: FileMapper.hh:110
General exception with possible name of object and message.
const char * what() const noexcept override
Override message in std::exception.
bool IsREBDSIMOrCombineOutputFile(TFile *file)
Definition: FileMapper.cc:116
void WarningMissingHistogram(const std::string &histName, const std::string &fileName)
Common print out method.
Definition: FileMapper.cc:147