21#include "BDSAcceleratorComponent.hh"
22#include "BDSBeamPipeInfo.hh"
23#include "BDSBendBuilder.hh"
24#include "BDSComponentFactory.hh"
26#include "BDSException.hh"
27#include "BDSFieldType.hh"
28#include "BDSGlobalConstants.hh"
29#include "BDSIntegratorSet.hh"
30#include "BDSIntegratorType.hh"
31#include "BDSIntegratorSetType.hh"
32#include "BDSMagnetStrength.hh"
34#include "BDSMagnet.hh"
35#include "BDSMagnetOuterInfo.hh"
36#include "BDSUtilities.hh"
38#include "parser/element.h"
39#include "parser/elementtype.h"
41#include "G4Transform3D.hh"
49 return !(finiteAngle || finiteField);
57 G4double incomingFaceAngle,
58 G4double outgoingFaceAngle,
59 G4bool buildFringeFields,
63 const G4String baseName = elementName;
66 const G4double arcLength = element->
l * CLHEP::m;
67 const G4double angle = (*st)[
"angle"];
68 G4double fint = element->
fint;
69 G4double fintx = element->
fintx;
70 G4double hgap = element->
hgap * CLHEP::m;
76 G4bool buildFringeIncoming = buildFringeFields;
77 G4bool buildFringeOutgoing = buildFringeFields;
80 BDSFieldType dipoleFieldType = finiteK1 ? BDSFieldType::dipolequadrupole : BDSFieldType::dipole;
82 if (buildFringeFields)
87 {buildFringeIncoming =
false;}
89 {buildFringeOutgoing =
false;}
100 {buildFringeIncoming =
true;}
102 {buildFringeOutgoing =
true;}
111 {buildFringeIncoming =
false;}
113 {buildFringeOutgoing =
false;}
120 if (prevElement->
type == ElementType::_SBEND)
122 buildFringeIncoming =
false;
124 {
throw BDSException( __METHOD_NAME__, prevElement->name +
" e2 clashes with " + elementName +
" e1");}
129 if (nextElement->
type == ElementType::_SBEND)
130 {buildFringeOutgoing =
false;}
136 buildFringeIncoming =
false;
137 buildFringeOutgoing =
false;
148 buildFringeIncoming =
false;
149 buildFringeOutgoing =
false;
178 BDSFieldType::dipole,
184 BDSComponentFactory::ScalingFieldOuter(element));
201 G4double centralArcLength = arcLength;
202 G4double centralAngle = angle;
203 G4double oneFringeAngle = 0;
205 {oneFringeAngle = (thinElementArcLength / arcLength) * angle;}
207 if (buildFringeIncoming)
209 centralArcLength -= thinElementArcLength;
210 centralAngle -= oneFringeAngle;
212 if (buildFringeOutgoing)
214 centralArcLength -= thinElementArcLength;
215 centralAngle -= oneFringeAngle;
219 G4double semiAngle = centralAngle / (G4double) nSBends;
220 G4double semiArcLength = centralArcLength / (G4double) nSBends;
223 (*semiStrength)[
"angle"] = semiAngle;
224 (*semiStrength)[
"length"] = semiArcLength;
226 G4double zExtentIn = 0;
227 G4double zExtentOut = 0;
228 G4bool fadeIn =
true;
229 G4bool fadeOut =
true;
237 if (incomingFaceAngle > 0)
238 {zExtentIn = 0.5*horizontalWidth*std::tan(incomingFaceAngle - 0.5*std::abs(semiAngle));}
239 else if (incomingFaceAngle < 0)
240 {zExtentIn = 0.5*horizontalWidth*std::tan(0.5*std::abs(semiAngle) + incomingFaceAngle);}
241 if (outgoingFaceAngle > 0)
242 {zExtentOut = 0.5*horizontalWidth*std::tan(outgoingFaceAngle - 0.5*std::abs(semiAngle));}
243 else if (outgoingFaceAngle < 0)
244 {zExtentOut = 0.5*horizontalWidth*std::tan(0.5*std::abs(semiAngle) + outgoingFaceAngle);}
247 if (std::abs(zExtentIn) < semiArcLength/4.0)
249 if (std::abs(zExtentOut) < semiArcLength/4.0)
253 G4String centralName = baseName +
"_even_ang";
259 0.5*semiAngle, 0.5*semiAngle, bpInfo,
270 BDSFieldType::dipole,
276 BDSComponentFactory::ScalingFieldOuter(element));
278 mgInfo->name = centralName;
294 G4double minimalRadius = 2*magnetOuterInfoCheck->MinimumIntersectionRadiusRequired();
296 if (!magnetOuterInfoCheck->hStyle)
298 switch (magnetOuterInfoCheck->geometryType.underlying())
300 case BDSMagnetGeometryType::polescircular:
301 case BDSMagnetGeometryType::polesfacet:
302 case BDSMagnetGeometryType::polesfacetcrop:
303 case BDSMagnetGeometryType::polessquare:
305 minimalRadius *= element->yokeOnInside ? 2.0 : 0.5;
313 delete magnetOuterInfoCheck;
316 if (buildFringeIncoming)
318 G4double e1 = element->
e1;
326 BDSMagnetStrength* fringeStIn = BDS::GetFringeMagnetStrength(element, st, oneFringeAngle,
327 e1, element->
e2, fintx,
true);
328 G4String segmentName = baseName +
"_e1_fringe";
329 G4double fringeAngleIn = 0.5*oneFringeAngle - incomingFaceAngle;
330 G4double fringeAngleOut = 0.5*oneFringeAngle + incomingFaceAngle;
332 segmentName, fringeStIn, brho,
333 integratorSet, dipoleFieldType);
343 G4double segmentAngleIn = 0;
344 G4double segmentAngleOut = 0;
345 G4int numberOfUniqueComponents = 1;
347 G4bool centralWedgeUsed =
false;
348 const G4int numSegmentsEitherSide = (nSBends - 1) / 2;
349 for (G4int i = 0; i < nSBends; ++i)
351 G4String name = baseName;
352 if (i < numSegmentsEitherSide)
355 {oneBend = centralWedge;}
358 name +=
"_"+std::to_string(numberOfUniqueComponents);
359 numberOfUniqueComponents++;
360 BDS::UpdateSegmentAngles(i,nSBends,semiAngle,incomingFaceAngle,outgoingFaceAngle,segmentAngleIn,segmentAngleOut);
362 segmentAngleIn, segmentAngleOut, semiStrength,
363 brho, integratorSet, yokeOnLeft, semiOuterField);
369 name +=
"_"+std::to_string(numberOfUniqueComponents);
370 numberOfUniqueComponents++;
371 segmentAngleIn = 0.5*semiAngle - incomingFaceAngle;
372 segmentAngleOut = 0.5*semiAngle;
374 segmentAngleIn, segmentAngleOut, semiStrength,
375 brho, integratorSet, yokeOnLeft, semiOuterField);
378 {oneBend = centralWedge;}
381 else if (i > numSegmentsEitherSide)
384 {oneBend = centralWedge;}
387 name +=
"_"+std::to_string(numberOfUniqueComponents);
388 numberOfUniqueComponents++;
389 BDS::UpdateSegmentAngles(i,nSBends,semiAngle,incomingFaceAngle,outgoingFaceAngle,segmentAngleIn,segmentAngleOut);
391 segmentAngleIn, segmentAngleOut, semiStrength,
392 brho, integratorSet, yokeOnLeft, semiOuterField);
396 if (i == (nSBends-1))
398 name +=
"_"+std::to_string(numberOfUniqueComponents);
399 numberOfUniqueComponents++;
400 segmentAngleIn = 0.5*semiAngle;
401 segmentAngleOut = 0.5*semiAngle - outgoingFaceAngle;
403 segmentAngleIn, segmentAngleOut, semiStrength,
404 brho, integratorSet, yokeOnLeft, semiOuterField);
407 {oneBend = centralWedge;}
411 {oneBend = centralWedge;}
416 centralWedgeUsed = centralWedgeUsed || (oneBend == centralWedge);
419 if (!centralWedgeUsed)
420 {
delete centralWedge;}
423 if (buildFringeOutgoing)
425 G4double e2 = element->
e2;
433 BDSMagnetStrength* fringeStOut = BDS::GetFringeMagnetStrength(element, st, oneFringeAngle,
434 element->
e1, e2, fintx,
false);
435 G4double fringeAngleIn = 0.5*oneFringeAngle + outgoingFaceAngle;
436 G4double fringeAngleOut = 0.5*oneFringeAngle - outgoingFaceAngle;
437 G4String segmentName = baseName +
"_e2_fringe";
440 segmentName, fringeStOut, brho,
441 integratorSet, dipoleFieldType);
447void BDS::UpdateSegmentAngles(G4int index,
450 G4double incomingFaceAngle,
451 G4double outgoingFaceAngle,
452 G4double& segmentAngleIn,
453 G4double& segmentAngleOut)
457 G4bool central = index == (nSBends - 1) / 2;
460 segmentAngleIn = 0.5*semiAngle;
461 segmentAngleOut = 0.5*semiAngle;
465 G4bool firstHalf = index + 1 < (nSBends + 1) / 2;
466 G4int numberToFadeOver = (nSBends - 1) / 2;
469 G4double delta = incomingFaceAngle / (G4double) numberToFadeOver;
470 G4double inputFaceAngle = incomingFaceAngle - ((G4double)index * delta);
471 G4double outputFaceAngle = incomingFaceAngle - ((G4double)(index+1) * delta);
473 segmentAngleIn = 0.5*semiAngle - inputFaceAngle;
474 segmentAngleOut = 0.5*semiAngle + outputFaceAngle;
478 G4double delta = outgoingFaceAngle / (G4double) numberToFadeOver;
479 G4int secondHalfIndex = index - ((nSBends + 1) / 2);
481 G4double inputFaceAngle = outgoingFaceAngle - ( (G4double)(numberToFadeOver - secondHalfIndex) * delta);
482 G4double outputFaceAngle = outgoingFaceAngle - ( (G4double)(numberToFadeOver - (secondHalfIndex + 1)) * delta);
484 segmentAngleIn = 0.5*semiAngle + inputFaceAngle;
485 segmentAngleOut = 0.5*semiAngle - outputFaceAngle;
490 const G4String& name,
506 magnetOuterInfo->name = name;
511 BDSFieldType dipoleFieldType = finiteK1 ? BDSFieldType::dipolequadrupole : BDSFieldType::dipole;
546 G4double incomingFaceAngle,
547 G4double outgoingFaceAngle,
548 G4bool buildFringeFields)
550 const G4String name = elementName;
555 G4bool buildFringeIncoming = buildFringeFields;
556 G4bool buildFringeOutgoing = buildFringeFields;
557 G4double fint = element->
fint;
558 G4double fintx = element->
fintx;
559 G4double hgap = element->
hgap * CLHEP::m;
562 const G4double angle = (*st)[
"angle"];
563 const G4double arcLength = (*st)[
"length"];
570 BDSFieldType dipoleFieldType = finiteK1 ? BDSFieldType::dipolequadrupole : BDSFieldType::dipole;
574 G4double incomingFaceAngleWRTSBend = incomingFaceAngle + angle*0.5;
575 G4double outgoingFaceangleWRTSBend = outgoingFaceAngle + angle*0.5;
577 {buildFringeIncoming =
false;}
579 {buildFringeOutgoing =
false;}
590 {buildFringeIncoming =
true;}
592 {buildFringeOutgoing =
true;}
595 G4double trackingPolefaceAngleIn = element->
e1;
596 G4double trackingPolefaceAngleOut = element->
e2;
599 trackingPolefaceAngleIn += angle*0.5;
600 trackingPolefaceAngleOut += angle*0.5;
603 G4double e1 = -incomingFaceAngle;
604 G4double e2 = -outgoingFaceAngle;
610 throw BDSException(__METHOD_NAME__, prevElement->name +
" has finite e2!\n Clashes with " + elementName +
" with finite e1");
612 if (prevElement->
type == ElementType::_RBEND)
614 buildFringeIncoming =
false;
622 throw BDSException(__METHOD_NAME__ + nextElement->name +
" has finite e1!\n Clashes with " + elementName +
" with finite e2");
624 if (nextElement->
type == ElementType::_RBEND)
626 buildFringeOutgoing =
false;
633 buildFringeIncoming =
false;
634 buildFringeOutgoing =
false;
648 G4double angleIn = incomingFaceAngle;
649 G4double angleOut = outgoingFaceAngle;
650 G4double fringeInOutputAngle = 0;
651 G4double centralInputFaceAngle = angleIn;
652 G4double centralOutputFaceAngle = angleOut;
653 G4double fringeOutInputAngle = 0;
654 G4double centralArcLength = arcLength;
655 G4double centralAngle = angle;
656 G4double oneFringeAngle = 0;
658 {oneFringeAngle = (thinElementArcLength / arcLength) * angle;}
660 if (buildFringeIncoming && buildFringeOutgoing)
662 centralArcLength -= 2*thinElementArcLength;
663 centralAngle -= 2*oneFringeAngle;
664 angleIn = e1 + (0.5*oneFringeAngle - 0.5*angle);
665 fringeInOutputAngle = -angleIn;
666 centralInputFaceAngle = e1;
667 centralOutputFaceAngle = e2;
668 fringeOutInputAngle = - (e2 + (0.5*oneFringeAngle - 0.5*angle));
669 angleOut = -fringeOutInputAngle;
671 else if (buildFringeIncoming)
673 centralArcLength -= thinElementArcLength;
674 centralAngle -= oneFringeAngle;
675 angleIn = e1 + (0.5*oneFringeAngle - 0.5*angle);
676 fringeInOutputAngle = -angleIn;
677 centralInputFaceAngle = e1 + 0.5*oneFringeAngle;
678 centralOutputFaceAngle = e2 - 0.5*oneFringeAngle;
680 else if (buildFringeOutgoing)
682 centralArcLength -= thinElementArcLength;
683 centralAngle -= oneFringeAngle;
684 centralInputFaceAngle = e1 - 0.5*oneFringeAngle;
685 centralOutputFaceAngle = e2 + 0.5*oneFringeAngle;
686 fringeOutInputAngle = - (e2 + (0.5*oneFringeAngle - 0.5*angle));;
687 angleOut = e2 + (0.5*oneFringeAngle - 0.5*angle);;
693 centralInputFaceAngle = e1;
694 centralOutputFaceAngle = e2;
697 if (buildFringeIncoming)
699 G4double trackingPolefaceAngle = trackingPolefaceAngleIn;
705 {trackingPolefaceAngle -= element->
e1;}
707 BDSMagnetStrength* fringeStIn = BDS::GetFringeMagnetStrength(element, st, oneFringeAngle,
708 trackingPolefaceAngle, trackingPolefaceAngleOut,
710 G4String fringeName = name +
"_e1_fringe";
716 integratorSet, dipoleFieldType);
721 (*st)[
"length"] = centralArcLength;
722 (*st)[
"angle"] = centralAngle;
729 mgInfo->name = elementName;
739 BDSFieldType::dipole,
745 BDSComponentFactory::ScalingFieldOuter(element));
752 elementName+
"_centre",
763 if (buildFringeOutgoing)
765 G4double trackingPolefaceAngle = trackingPolefaceAngleOut;
771 {trackingPolefaceAngle -= element->
e2;}
772 BDSMagnetStrength* fringeStOut = BDS::GetFringeMagnetStrength(element, st, oneFringeAngle,
773 trackingPolefaceAngleIn, trackingPolefaceAngle,
775 G4String fringeName = name +
"_e2_fringe";
780 integratorSet, dipoleFieldType);
790 const G4String& name,
797 beamPipeInfo->
beamPipeType = BDSBeamPipeType::circularvacuum;
799 angleIn, angleOut, beamPipeInfo);
800 magnetOuterInfo->geometryType = BDSMagnetGeometryType::none;
801 magnetOuterInfo->name = name;
802 magnetOuterInfo->buildEndPieces =
false;
814 return new BDSMagnet(BDSMagnetType::dipolefringe,
827 G4double incomingFaceAngle,
828 G4double outgoingFaceAngle,
829 G4double aperturePrecision)
836 G4double totalAngle = std::abs(angle) + std::abs(incomingFaceAngle) + std::abs(outgoingFaceAngle);
837 G4int nSBends = (G4int) ceil(std::sqrt(length*totalAngle/2/aperturePrecision));
842 if (nSBends % 2 == 0)
845 G4cout << __METHOD_NAME__ <<
" splitting sbend into " << nSBends <<
" sbends" << G4endl;
858 {intType = integratorSet->
Integrator(BDSFieldType::dipolequadrupole);}
868 G4double fringeAngle,
876 (*fringeSt)[
"angle"] = fringeAngle;
877 (*fringeSt)[
"e1"] = e1;
878 (*fringeSt)[
"e2"] = e2;
879 (*fringeSt)[
"fint"] = element->
fint;
880 (*fringeSt)[
"fintx"] = fintx;
881 (*fringeSt)[
"fintk2"] = element->
fintK2;
882 (*fringeSt)[
"fintk2"] = element->
fintxK2;
883 (*fringeSt)[
"hgap"] = element->
hgap * CLHEP::m;
884 (*fringeSt)[
"isentrance"] = isEntrance;
885 (*fringeSt)[
"h1"] = element->
h1 / CLHEP::m;
886 (*fringeSt)[
"h2"] = element->
h2 / CLHEP::m;
Abstract class that represents a component of an accelerator.
Holder class for all information required to describe a beam pipe model.
BDSBeamPipeType beamPipeType
Public member for direct access.
static G4bool YokeOnLeft(const GMAD::Element *el, const BDSMagnetStrength *st)
static G4double PrepareHorizontalWidth(GMAD::Element const *el, G4double defaultHorizontalWidth=-1)
Prepare the element horizontal width in Geant4 units - if not set, use the global default.
static G4Transform3D CreateFieldTransform(const BDSTiltOffset *tiltOffset)
Create a transform from a tilt offset. If nullptr, returns identity transform.
static BDSMagnetOuterInfo * PrepareMagnetOuterInfo(const G4String &elementNameIn, const GMAD::Element *el, const BDSMagnetStrength *st, const BDSBeamPipeInfo *beamPipe, G4double defaultHorizontalWidth=-1, G4double defaultVHRatio=1.0, G4double defaultCoilWidthFraction=-1, G4double defaultCoilHeightFraction=-1)
static void CheckBendLengthAngleWidthCombo(G4double arcLength, G4double angle, G4double horizontalWidth, const G4String &name="not given")
static BDSFieldInfo * PrepareMagnetOuterFieldInfo(const BDSMagnetStrength *vacuumSt, const BDSFieldType &fieldType, const BDSBeamPipeInfo *bpInfo, const BDSMagnetOuterInfo *outerInfo, const G4Transform3D &fieldTransform, const BDSIntegratorSet *integratorSetIn, G4double brhoIn, G4double outerFieldScaling=1.0)
Prepare the field definition for the yoke of a magnet.
static BDSBeamPipeInfo * PrepareBeamPipeInfo(GMAD::Element const *el, const G4ThreeVector &inputFaceNormal=G4ThreeVector(0, 0,-1), const G4ThreeVector &outputFaceNormal=G4ThreeVector(0, 0, 1))
General exception with possible name of object and message.
All info required to build complete field of any type.
static BDSGlobalConstants * Instance()
Access method.
Which integrator to use for each type of magnet / field object.
G4bool IsMatrixIntegratorSet() const
Accessor for bool of is the integrator set matrix style.
BDSIntegratorType Integrator(const BDSFieldType field) const
Get appropriate integrator based on the field type.
A class that hold multiple accelerator components.
void AddComponent(BDSAcceleratorComponent *component)
Add a component to the line.
Efficient storage of magnet strengths.
Abstract base class that implements features common to all magnets.
G4bool ZeroStrengthDipole(const BDSMagnetStrength *st)
Return whether finite angle or field for a dipole.
G4int CalculateNSBendSegments(G4double length, G4double angle, G4double incomingFaceAngle=0, G4double outgoingFaceAngle=0, G4double aperturePrecision=1.0)
BDSMagnet * BuildDipoleFringe(const GMAD::Element *element, G4double angleIn, G4double angleOut, const G4String &name, BDSMagnetStrength *st, G4double brho, const BDSIntegratorSet *integratorSet, BDSFieldType dipoleFieldType)
BDSLine * BuildRBendLine(const G4String &elementName, const GMAD::Element *element, const GMAD::Element *prevElement, const GMAD::Element *nextElement, G4double brho, BDSMagnetStrength *st, const BDSIntegratorSet *integratorSet, G4double incomingFaceAngle, G4double outgoingFaceAngle, G4bool buildFringeFields)
BDSAcceleratorComponent * BuildSBendLine(const G4String &elementName, const GMAD::Element *element, BDSMagnetStrength *st, G4double brho, const BDSIntegratorSet *integratorSet, G4double incomingFaceAngle, G4double outgoingFaceAngle, G4bool buildFringeFields, const GMAD::Element *prevElement, const GMAD::Element *nextElement)
BDSIntegratorType GetDipoleIntegratorType(const BDSIntegratorSet *integratorSet, const GMAD::Element *element)
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())
BDSMagnet * BuildSingleSBend(const GMAD::Element *element, const G4String &name, G4double arcLength, G4double angle, G4double angleIn, G4double angleOut, const BDSMagnetStrength *strength, G4double brho, const BDSIntegratorSet *integratorSet, G4bool yokeOnLeft, const BDSFieldInfo *outerFieldIn)
Function to return a single sector bend section.
Parser namespace for GMAD language. Combination of Geant4 and MAD.
double hgap
half distance of pole separation for purposes of fringe fields - 'half gap'
double fintx
fringe field integral at the dipole exit
double h2
output pole face curvature for bends
double h1
input pole face curvature for bends
double fintxK2
second fringe field integral at the dipole exit - for TRANSPORT matching
double fint
fringe field integral at the dipole entrance
double e2
output pole face rotation for bends
double e1
input pole face rotation for bends
double fintK2
second fringe field integral at the dipole entrance - for TRANSPORT matching
ElementType type
element enum