BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
bdsimCombine.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
25#include "BDSOutputROOTEventHeader.hh"
26
27#include "TChain.h"
28#include "TFile.h"
29#include "TTree.h"
30
31#include <exception>
32#include <glob.h>
33#include <iostream>
34#include <string>
35#include <vector>
36
37int main(int argc, char* argv[])
38{
39 if (argc < 3)
40 {
41 std::cout << "usage: bdsimCombine result.root file1.root file2.root ..." << std::endl;
42 exit(1);
43 }
44
45 // build input file list
46 std::vector<std::string> inputFiles;
47 for (int i = 2; i < argc; ++i)
48 {inputFiles.emplace_back(std::string(argv[i]));}
49 // see if we're globbing files
50 if (inputFiles[0].find('*') != std::string::npos)
51 {
52 std::vector<std::string> fileNamesTemp;
53 glob_t glob_result;
54 glob(inputFiles[0].c_str(),GLOB_TILDE,nullptr,&glob_result);
55 for (unsigned int i=0; i<glob_result.gl_pathc; ++i)
56 {fileNamesTemp.emplace_back(glob_result.gl_pathv[i]);}
57 globfree(&glob_result);
58 inputFiles = fileNamesTemp;
59 }
60
61 // checks
62 if (inputFiles.size() == 1)
63 {
64 std::cout << "Only one input file provided \"" << inputFiles[0] << "\" - no point." << std::endl;
65 exit(1);
66 }
67 // check for wrong order of arguments which is common mistake
68 std::string outputFile = std::string(argv[1]);
69 if (outputFile.find('*') != std::string::npos)
70 {
71 std::cerr << "First argument for output file \"" << outputFile << "\" contains an *." << std::endl;
72 std::cerr << "Should only be a singular file - check order of arguments." << std::endl;
73 exit(1);
74 }
75
76 // merge event trees - this is done by ROOT and it writes the output file and closes it
77 TChain* eventsMerged = new TChain("Event");
78 for (const auto& filename : inputFiles)
79 {eventsMerged->Add(filename.c_str());}
80 std::cout << "Beginning merge of Event Tree" << std::endl;
81 Long64_t operationCode = eventsMerged->Merge(outputFile.c_str());
82 if (operationCode == 0)
83 {// from inspection of ROOT TChain.cxx ~line 1866, it returns 0 if there's a problem
84 std::cerr << "Problem in TTree::Merge opening output file \"" << outputFile << "\"" << std::endl;
85 return 1;
86 }
87 else
88 {std::cout << "Finished merge of Event Tree" << std::endl;}
89
90 // loop over input files loading headers to accumulate number of original events
91 unsigned long long int nOriginalEvents = 0;
92 bool skimmedFile = false;
93 unsigned long int i = 0;
94 std::cout << "Counting number of original events from headers of files" << std::endl;
95 for (const auto& filename : inputFiles)
96 {
97 TFile* f = new TFile(filename.c_str(), "READ");
99 {
100 std::cout << "File \"" << filename << "\" skipped as not a valid BDSIM file" << std::endl;
101 if (!f->IsZombie())
102 {f->Close();}
103 delete f;
104 continue;
105 }
106 std::cout << "Accumulating> " << filename << std::endl;
107 TTree* headerTree = dynamic_cast<TTree*>(f->Get("Header")); // should be safe given check we've just done
108 if (!headerTree)
109 {
110 std::cerr << "Problem getting header from file " << filename << std::endl;
111 f->Close();
112 delete f;
113 continue;
114 }
115 Header* headerLocal = new Header();
116 headerLocal->SetBranchAddress(headerTree);
117 headerTree->GetEntry(0);
118 skimmedFile = skimmedFile || headerLocal->header->skimmedFile;
119 if (headerLocal->header->skimmedFile)
120 {nOriginalEvents += headerLocal->header->nOriginalEvents;}
121 else
122 {// unskimmed file which won't record the number of events in the header, so we inspect the Event Tree
123 TTree* eventTree = dynamic_cast<TTree*>(f->Get("Event"));
124 if (eventTree)
125 {
126 Long64_t nEntries = eventTree->GetEntries();
127 nOriginalEvents += (unsigned long long int)nEntries;
128 }
129 else
130 {std::cerr << "Problem getting Event tree in file " << filename << std::endl;}
131 }
132 delete headerLocal;
133 f->Close();
134 delete f;
135 i++;
136 }
137
138 // checks
139 if (i == 0)
140 {std::cerr << "No valid files found" << std::endl; return 1;}
141
142 // now we produce a new header and update the file as well as copy over the other trees from the first valid
143 // input file in the list (i.e. tolerate the odd zombie file from a big run)
144 TFile* input = nullptr;
145 bool validFile = false;
146 i = 0;
147 while (!validFile)
148 {// assume we have one valid file from check above
149 if (i > (unsigned long int)inputFiles.size())
150 {return 1;}
151 input = new TFile(inputFiles[i].c_str(), "READ");
152 validFile = RBDS::IsBDSIMOutputFile(input);
153 i++;
154 }
155
156 TFile* output = new TFile(outputFile.c_str(), "UPDATE");
157 if (output->IsZombie())
158 {std::cerr << "Could not reopen output file to add other trees"; return 1;}
159 output->cd();
161 headerOut->Fill(std::vector<std::string>(), inputFiles); // updates time stamp
162 headerOut->SetFileType("BDSIM");
163 headerOut->skimmedFile = skimmedFile;
164 headerOut->nOriginalEvents = nOriginalEvents;
165 TTree* headerTree = new TTree("Header", "BDSIM Header");
166 headerTree->Branch("Header.", "BDSOutputROOTEventHeader", headerOut);
167 headerTree->Fill();
168 output->Write(nullptr,TObject::kOverwrite);
169
170 // go over all other trees and copy them (in the original order) from the first file to the output
171 std::cout << "Merging rest of file contents" << std::endl;
172 std::vector<std::string> treeNames = {"ParticleData", "Beam", "Options", "Model", "Run"};
173 for (const auto& tn : treeNames)
174 {
175 TTree* original = dynamic_cast<TTree*>(input->Get(tn.c_str()));
176 if (!original)
177 {
178 std::cerr << "Failed to load Tree named " << tn << std::endl;
179 delete output;
180 delete input;
181 return 1;
182 }
183 auto clone = original->CloneTree();
184 clone->AutoSave();
185 }
186
187 output->Close();
188 delete output;
189
190 std::cout << "Combined result of " << inputFiles.size() << " files written to: " << outputFile << std::endl;
191 std::cout << "Run histograms are not summed" << std::endl; // TODO
192 return 0;
193}
Information about the software and the file.
void SetFileType(std::string fileTypeIn)
Update the file type.
bool skimmedFile
Whether this is a skimmed output file.
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
bool IsBDSIMOutputFile(TFile *file, int *dataVersion=nullptr)
Definition: FileMapper.cc:66