BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
rebdsimOrbit.cc
Go to the documentation of this file.
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2023.
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*/
23#include <iostream>
24#include <string>
25
26#include "DataLoader.hh"
27#include "EventAnalysisOrbit.hh"
28#include "RBDSException.hh"
29
30#include "BDSOutputROOTEventHeader.hh"
31
32#include "TFile.h"
33#include "TTree.h"
34
35void usage()
36{
37 std::cout << "usage: rebdsimOrbit <dataFile> <outputfile> <index>" << std::endl;
38 std::cout << " <datafile> - root file to operate on ie run1.root" << std::endl;
39 std::cout << " <outputfile> - name of output file ie optics.root" << std::endl;
40 std::cout << " <index> - index of orbit to pull from file (default 1)" << std::endl;
41 std::cout << "This only works on a single file - no wildcards." << std::endl;
42}
43
44int main(int argc, char* argv[])
45{
46 if (argc < 3 || argc >= 5)
47 {
48 usage();
49 return 1;
50 }
51
52 std::string inputFileName = std::string(argv[1]);
53 std::string outputFileName = std::string(argv[2]);
54 int index = 0;
55 try
56 {
57 if (argc == 4)
58 {index = std::stoi(argv[3]);}
59
60 DataLoader dl = DataLoader(inputFileName, false, true, true);
61 EventAnalysisOrbit* evtAnalysis = new EventAnalysisOrbit(dl.GetEvent(),
62 dl.GetEventTree(),
63 true, true);
64 evtAnalysis->ExtractOrbit(index);
65
66 TFile* outputFile = new TFile(outputFileName.c_str(), "RECREATE");
67
68 // add header for file type and version details
69 outputFile->cd();
71 headerOut->Fill(dl.GetFileNames()); // updates time stamp
72 headerOut->SetFileType("REBDSIM");
73 TTree* headerTree = new TTree("Header", "REBDSIM Header");
74 headerTree->Branch("Header.", "BDSOutputROOTEventHeader", headerOut);
75 headerTree->Fill();
76 headerTree->Write();
77
78 evtAnalysis->WriteOrbit(outputFile);
79 outputFile->Close();
80 std::cout << "Result written to: " << outputFileName << std::endl;
81 delete outputFile;
82 delete evtAnalysis;
83 }
84 catch (const RBDSException& error)
85 {std::cerr << error.what() << std::endl; return 1;}
86 catch (const std::exception& error)
87 {std::cerr << error.what() << std::endl; return 1;}
88
89 return 0;
90}
Information about the software and the file.
void SetFileType(const std::string &fileTypeIn)
Update the file type.
void Fill(const std::vector< std::string > &analysedFilesIn=std::vector< std::string >(), const std::vector< std::string > &combinedFilesIn=std::vector< std::string >())
Loader for a ROOT file using classes used to generate the file.
Definition: DataLoader.hh:46
std::vector< std::string > GetFileNames()
Accessor.
Definition: DataLoader.hh:81
Event * GetEvent()
Accessor.
Definition: DataLoader.hh:91
TChain * GetEventTree()
Accessor.
Definition: DataLoader.hh:98
Simple analysis to pull out first hit in each sampler.
void ExtractOrbit(int index)
Extract an orbit from the data.
void WriteOrbit(TFile *f)
Write orbit to a ROOT file.
General exception with possible name of object and message.
const char * what() const noexcept override
Override message in std::exception.