19#include "BDSOutputROOTEventModel.hh"
22#include "BDSAcceleratorModel.hh"
23#include "BDSBeamline.hh"
24#include "BDSBeamlineElement.hh"
25#include "BDSBeamPipeInfo.hh"
26#include "BDSMagnet.hh"
27#include "BDSMagnetStrength.hh"
28#include "BDSPhysicalVolumeInfoRegistry.hh"
29#include "BDSSamplerRegistry.hh"
31#include "G4RotationMatrix.hh"
49 storeCollimatorInfo(false),
50 storeCavityInfo(false),
66 std::vector<std::pair<float, int> > distanceAndIndex;
67 distanceAndIndex.reserve(midRefPos.size());
68 for (
int i = 0; i < (int)midRefPos.size(); i++)
69 {distanceAndIndex.emplace_back(std::make_pair((midRefPos[i]-point).Mag(), i));}
72 {
bool operator()(std::pair<float, int>& a, std::pair<float, int>& b)
const {
return a.first < b.first;};};
73 std::sort(distanceAndIndex.begin(), distanceAndIndex.end(), customLess());
74 return distanceAndIndex[0].second;
80 componentName.clear();
81 placementName.clear();
82 componentType.clear();
103 beamPipeType.clear();
104 beamPipeAper1.clear();
105 beamPipeAper2.clear();
106 beamPipeAper3.clear();
107 beamPipeAper4.clear();
145 pvNameWPointer.clear();
159 scoringMeshTranslation.clear();
160 scoringMeshRotation.clear();
161 scoringMeshName.clear();
163 materialIDToName.clear();
164 materialNameToID.clear();
166 samplerNamesUnique.clear();
167 samplerSPosition.clear();
168 samplerCNamesUnique.clear();
169 samplerSNamesUnique.clear();
170 samplerCRadius.clear();
171 samplerSRadius.clear();
176 G4bool storeCavityInfoIn):
178 storeCollimatorInfo(storeCollimatorInfoIn),
179 storeCavityInfo(storeCavityInfoIn),
186 TRotation rr = TRotation();
192 CLHEP::Hep3Vector axis;
194 rm->getAngleAxis(rotAngle, axis);
195 rr.Rotate(rotAngle, TVector3(axis.x(), axis.y(), axis.z()));
206 return TVector3(v.x() / CLHEP::m, v.y() / CLHEP::m, v.z() / CLHEP::m);
210 const std::map<G4String, G4int>& collimatorIndicesByNameIn,
211 const std::vector<BDSOutputROOTEventCollimatorInfo>& collimatorInfoIn,
212 const std::vector<G4String>& collimatorBranchNamesIn,
213 const std::vector<G4int>& cavityIndicesIn,
214 const std::map<G4String, G4int>& cavityIndicesByNameIn,
215 const std::vector<BDSOutputROOTEventCavityInfo>& cavityInfoIn,
216 const std::vector<G4String>& cavityBranchNamesIn,
217 const std::map<G4String, G4Transform3D>* scorerMeshPlacements,
218 const std::map<short int, G4String>* materialIDToNameUnique,
219 G4bool storeTrajectory)
222 for (
const auto& nameSPos : sr->GetUniquePlaneNamesAndSPosition())
224 samplerNamesUnique.push_back(std::string(nameSPos.first) +
".");
225 samplerSPosition.push_back((
double) nameSPos.second / CLHEP::m);
227 for (
const auto& name : sr->GetUniqueNamesCylinder())
228 {samplerCNamesUnique.push_back(std::string(name) +
".");}
229 for (
const auto& name : sr->GetUniqueNamesSphere())
230 {samplerSNamesUnique.push_back(std::string(name) +
".");}
231 samplerCRadius = sr->GetUniqueNameToRadiusCylinder();
232 samplerSRadius = sr->GetUniqueNameToRadiusSphere();
234 for (
const auto& name : collimatorBranchNamesIn)
236 for (
const auto& name : cavityBranchNamesIn)
239 if (scorerMeshPlacements)
241 for (
const auto& kv : *scorerMeshPlacements)
243 const G4String& name = kv.first;
244 scoringMeshTranslation[name] =
ConvertToROOT(kv.second.getTranslation());
245 scoringMeshRotation[name] =
ConvertToROOT(kv.second.getRotation());
246 scoringMeshName.emplace_back(name);
256 for (
const auto value : collimatorIndicesIn)
261 for (
const auto& kv : collimatorIndicesByNameIn)
269 for (
const auto value : cavityIndicesIn)
274 for (
const auto& kv : cavityIndicesByNameIn)
280 if (materialIDToNameUnique && storeTrajectory)
282 for (
const auto& kv : *materialIDToNameUnique)
284 materialIDToName[kv.first] = (std::string)kv.second;
285 materialNameToID[(std::string)kv.second] = kv.first;
289 n = (int)beamline->
size();
291 for (
auto i = beamline->
begin(); i != beamline->
end(); ++i)
293 componentName.push_back((*i)->GetName());
294 placementName.push_back((*i)->GetPlacementName());
295 componentType.push_back((*i)->GetType());
296 length.push_back((
float)((*i)->GetAcceleratorComponent()->GetArcLength() / CLHEP::m));
297 angle.push_back((
float)((*i)->GetAcceleratorComponent()->GetAngle() / CLHEP::radian));
304 staRefPos.emplace_back(
ConvertToROOT((*i)->GetReferencePositionStart()));
305 midRefPos.emplace_back(
ConvertToROOT((*i)->GetReferencePositionMiddle()));
306 endRefPos.emplace_back(
ConvertToROOT((*i)->GetReferencePositionEnd()));
307 staRefRot.push_back(
ConvertToROOT((*i)->GetReferenceRotationStart()));
308 midRefRot.push_back(
ConvertToROOT((*i)->GetReferenceRotationMiddle()));
309 endRefRot.push_back(
ConvertToROOT((*i)->GetReferenceRotationEnd()));
315 tilt.push_back((
float)(to->
GetTilt() / CLHEP::rad));
316 offsetX.push_back((
float)(to->
GetXOffset() / CLHEP::m));
317 offsetY.push_back((
float)(to->
GetYOffset() / CLHEP::m));
322 offsetX.push_back(0);
323 offsetY.push_back(0);
327 staS.push_back((
float)((*i)->GetSPositionStart() / CLHEP::m));
328 midS.push_back((
float)((*i)->GetSPositionMiddle() / CLHEP::m));
329 endS.push_back((
float)((*i)->GetSPositionEnd() / CLHEP::m));
333 beamPipeType.push_back(beampipeinfo ? beampipeinfo->
beamPipeType.ToString() :
"");
334 beamPipeAper1.push_back(beampipeinfo ? beampipeinfo->
aper1 / CLHEP::m : 0);
335 beamPipeAper2.push_back(beampipeinfo ? beampipeinfo->
aper2 / CLHEP::m : 0);
336 beamPipeAper3.push_back(beampipeinfo ? beampipeinfo->
aper3 / CLHEP::m : 0);
337 beamPipeAper4.push_back(beampipeinfo ? beampipeinfo->
aper4 / CLHEP::m : 0);
340 const auto accComp = (*i)->GetAcceleratorComponent();
341 material.push_back(accComp->Material());
344 std::vector<std::vector<float>*> localNorm = {&k1,&k2,&k3,&k4,&k5,&k6,&k7,&k8,&k9,&k10,&k11,&k12};
345 std::vector<std::vector<float>*> localSkew = {&k1s,&k2s,&k3s,&k4s,&k5s,&k6s,&k7s,&k8s,&k9s,&k10s,&k11s,&k12s};
350 for (
int j = 0; j < (int)localNorm.size(); j++)
351 {localNorm[j]->push_back(0);}
352 for (
int j = 0; j < (int)localSkew.size(); j++)
353 {localSkew[j]->push_back(0);}
365 fintxk2.push_back(0);
370 std::vector<std::string> localPVNames;
371 std::vector<std::string> localPVNamesWPointer;
374 auto name = [](G4VPhysicalVolume* pv){
return pv->GetName();};
375 std::transform(setOfPVs->begin(), setOfPVs->end(), std::back_inserter(localPVNames), name);
376 auto fullName = [](G4VPhysicalVolume* pv)
378 std::stringstream ss;
379 ss << static_cast<const void*>(pv);
380 std::string addressName = ss.str();
381 return pv->GetName() + addressName;
383 std::transform(setOfPVs->begin(), setOfPVs->end(), std::back_inserter(localPVNamesWPointer), fullName);
386 pvName.push_back(localPVNames);
387 pvNameWPointer.push_back(localPVNamesWPointer);
401 for (
int j = 0; j < (int)localNorm.size(); j++)
402 {localNorm[j]->push_back((
float)normComponents[j]);}
404 for (
int j = 0; j < (int)localSkew.size(); j++)
405 {localSkew[j]->push_back((
float)skewComponents[j]);}
415 fint.push_back((
float)(*ms)[
"fint"]);
416 fintx.push_back((
float)(*ms)[
"fintx"]);
417 fintk2.push_back((
float)(*ms)[
"fintk2"]);
418 fintxk2.push_back((
float)(*ms)[
"fintxk2"]);
const BDSBeamline * BeamlineMain() const
Accessor.
Holder class for all information required to describe a beam pipe model.
G4double aper3
Public member for direct access.
G4double aper1
Public member for direct access.
G4double aper4
Public member for direct access.
G4double aper2
Public member for direct access.
BDSBeamPipeType beamPipeType
Public member for direct access.
A vector of BDSBeamlineElement instances - a beamline.
BeamlineVector::size_type size() const
Get the number of elements.
iterator begin()
Iterator mechanics.
iterator end()
Iterator mechanics.
Efficient storage of magnet strengths.
std::vector< G4double > NormalComponents() const
Accessor for all normal components - k1 - k12.
std::vector< G4double > SkewComponents() const
Accessor for all skew components - k1 - k12.
static G4double Unit(const G4String &key)
Access a unit factor for a given key.
Abstract base class that implements features common to all magnets.
Information stored per model representing accelerator.
std::vector< BDSOutputROOTEventCollimatorInfo > collimatorInfo
Collimator information explicitly.
bool storeCollimatorInfo
Whether optional collimator information was stored.
std::vector< float > vkick
Vertical fractional momentum kick.
std::vector< float > eField
E field in V/m.
std::vector< std::string > collimatorBranchNamesUnique
Vector of all collimator branch names in event tree used to load data.
int nCollimators
Number of collimators in beam line.
int nCavities
Number of cavities in beam line.
virtual ~BDSOutputROOTEventModel()
Destructor.
int findNearestElement(const TVector3 &point) const
Find element index closest to point.
std::vector< float > bField
B field in T.
std::vector< float > hkick
Horizontal fractional momentum kick.
std::vector< int > collimatorIndices
std::vector< std::string > material
Material associated with element if any.
std::vector< BDSOutputROOTEventCavityInfo > cavityInfo
cavity information explicitly.
void Flush()
Initialise all members.
virtual void Fill(const std::vector< G4int > &collimatorIndicesIn={}, const std::map< G4String, G4int > &collimatorIndicesByNameIn={}, const std::vector< BDSOutputROOTEventCollimatorInfo > &collimatorInfoIn={}, const std::vector< G4String > &collimatorBranchNamesIn={}, const std::vector< G4int > &cavityIndicesIn={}, const std::map< G4String, G4int > &cavityIndicesByNameIn={}, const std::vector< BDSOutputROOTEventCavityInfo > &cavityInfoIn={}, const std::vector< G4String > &cavityBranchNamesIn={}, const std::map< G4String, G4Transform3D > *scorerMeshPlacements=nullptr, const std::map< short int, G4String > *materialIDToNameUnique=nullptr, G4bool storeTrajectory=false)
Fill root output.
std::vector< int > cavityIndices
std::vector< std::string > cavityBranchNamesUnique
Vector of all cavity branch names in event tree used to load data.
BDSOutputROOTEventModel()
Default constructor.
std::map< std::string, int > collimatorIndicesByName
Similar cache but by name of collimator as built by BDSIM.
std::vector< float > ks
Solenoid strength.
bool storeCavityInfo
Whether optional cavity information was stored.
TRotation ConvertToROOT(const G4RotationMatrix *rm) const
Utility function.
std::map< std::string, int > cavityIndicesByName
Similar cache but by name of cavity as built by BDSIM.
const std::set< G4VPhysicalVolume * > * PVsForBeamlineElement(BDSBeamlineElement *element) const
Access a set of volumes registered for the placement of a beamline element.
static BDSPhysicalVolumeInfoRegistry * Instance()
Singleton accessor.
static BDSSamplerRegistry * Instance()
Accessor for registry.
A holder for any placement offsets and rotations for a BDSAcceleratorComponent.
G4double GetYOffset() const
Accessor.
G4double GetXOffset() const
Accessor.
G4double GetTilt() const
Accessor.