19#include "BDSColourScaleMagma.hh"
20#include "BDSColourScaleViridis.hh"
21#include "BDSFieldQueryForVis.hh"
22#include "BDSFieldQueryInfo.hh"
23#include "BDSVisFieldModel.hh"
25#include "G4ArrowModel.hh"
30#include "G4ThreeVector.hh"
32#include "G4VGraphicsScene.hh"
33#include "G4VMarker.hh"
34#include "G4VisAttributes.hh"
35#include "BDSWarning.hh"
51BDSVisFieldModel::BDSVisFieldModel(
const std::vector<BDSFieldQueryInfo*>& queriesIn):
56 fGlobalTag =
"BDSVisFieldModel";
57 fGlobalDescription =
"Field view #" + std::to_string(instanceCounter);
61 pointVisB =
new G4VisAttributes(bFieldColour.GetValue(0));
63 pointVisE =
new G4VisAttributes(eFieldColour.GetValue(0));
65 pointVisB->SetForceSolid(
true);
66 pointVisE->SetForceSolid(
true);
69BDSVisFieldModel::~BDSVisFieldModel()
75void BDSVisFieldModel::DescribeYourselfTo(G4VGraphicsScene& sceneHandler)
81 for (
const auto& query :
queries)
84 G4double maxFieldB = querier.MaxFieldB();
85 G4double maxFieldE = querier.MaxFieldE();
86 const std::vector<std::array<G4double, 9>> xyzBEs = querier.Values();
89 if (arrowLength > 0.01 * std::numeric_limits<G4double>::max())
90 {arrowLength = 1e-3 * sceneHandler.GetExtent().GetExtentRadius();}
92 G4double arrowWidth = 0.3*arrowLength;
93 G4double pointSize = std::max(0.1*arrowWidth, 5.0);
95 G4bool drawArrows = query->drawArrows;
96 G4bool drawBoxes = query->drawBoxes;
98 if (!drawArrows && !drawBoxes)
100 BDS::Warning(
"Neither drawArrows or drawBoxes are turned on in query \"" + query->name +
"\" -> skipping");
104 G4double boxAlpha = query->boxAlpha;
107 BDS::Warning(
"boxAlpha in query \"" + query->name +
"\" must be > 0. Setting to 0.01");
112 BDS::Warning(
"boxAlpha in query \"" + query->name +
"\" must be <= 1. Setting to 1.0");
117 G4AffineTransform qTransform = query->globalTransform;
118 G4RotationMatrix qRotation = qTransform.NetRotation().inverse();
119 G4Box aBox(
"b", boxHalfSize.x(), boxHalfSize.y(), boxHalfSize.z());
120 G4Polyhedron* baseBoxPoly = aBox.GetPolyhedron();
122 if (query->queryMagnetic)
124 G4Colour boxZeroColourB =
pointVisB->GetColor();
125 G4Colour boxZeroColourB2(boxZeroColourB.GetRed(), boxZeroColourB.GetGreen(), boxZeroColourB.GetBlue(), boxAlpha);
126 G4VisAttributes* boxZeroVisB =
new G4VisAttributes(boxZeroColourB2);
127 G4String arrowPrefix = query->name +
"_B_";
128 for (
const auto& xyzBE: xyzBEs)
131 G4ThreeVector B(xyzBE[3], xyzBE[4], xyzBE[5]);
132 G4ThreeVector unitB = B.unit();
133 G4double bMag = B.mag();
134 G4ThreeVector midPoint(xyzBE[0], xyzBE[1], xyzBE[2]);
140 G4Polyhedron boxPoly(*baseBoxPoly);
141 boxPoly.Transform(G4Transform3D(qRotation, midPoint));
145 sceneHandler.BeginPrimitives();
148 G4Circle aPoint({xyzBE[0], xyzBE[1], xyzBE[2]});
150 aPoint.SetDiameter(G4VMarker::SizeType::world, pointSize);
151 aPoint.SetFillStyle(G4VMarker::FillStyle::filled);
152 sceneHandler.AddPrimitive(aPoint);
156 boxPoly.SetVisAttributes(boxZeroVisB);
157 sceneHandler.AddPrimitive(boxPoly);
159 sceneHandler.EndPrimitives();
162 G4ThreeVector startPoint = midPoint - 0.5*arrowLength*unitB;
163 G4ThreeVector endPoint = midPoint + 0.5*arrowLength*unitB;
164 G4double normalisedValue = bMag / maxFieldB;
165 if (!std::isfinite(normalisedValue))
166 {normalisedValue = 0.0;}
167 G4Colour arrowColour = bFieldColour.GetValue(normalisedValue);
170 G4Colour boxColour(arrowColour.GetRed(), arrowColour.GetGreen(), arrowColour.GetBlue(), boxAlpha);
171 sceneHandler.BeginPrimitives();
172 boxPoly.SetVisAttributes(G4VisAttributes(boxColour));
173 sceneHandler.AddPrimitive(boxPoly);
174 sceneHandler.EndPrimitives();
178 std::stringstream ss;
179 ss << std::setw(2) << std::fixed << xyzBE[0] <<
"_" << xyzBE[1] <<
"_" << xyzBE[2];
180 std::string arrowName = ss.str();
181 G4ArrowModel FArrow(startPoint.x(), startPoint.y(), startPoint.z(),
182 endPoint.x(), endPoint.y(), endPoint.z(),
183 arrowWidth, arrowColour, arrowName);
184 FArrow.DescribeYourselfTo(sceneHandler);
190 if (query->queryElectric)
192 G4Colour boxZeroColourE =
pointVisE->GetColor();
193 G4Colour boxZeroColourE2(boxZeroColourE.GetRed(),
194 boxZeroColourE.GetGreen(),
195 boxZeroColourE.GetBlue(),
197 G4VisAttributes* boxZeroVisE =
new G4VisAttributes(boxZeroColourE2);
198 G4String arrowPrefix = query->name +
"_E_";
199 for (
const auto& xyzBE: xyzBEs)
202 G4ThreeVector E(xyzBE[6], xyzBE[7], xyzBE[8]);
203 G4ThreeVector unitE = E.unit();
204 G4double eMag = E.mag();
205 G4ThreeVector midPoint(xyzBE[0], xyzBE[1], xyzBE[2]);
207 G4Polyhedron boxPoly(*baseBoxPoly);
208 boxPoly.Transform(G4Transform3D(qRotation, midPoint));
211 sceneHandler.BeginPrimitives();
214 G4Circle aPoint({xyzBE[0], xyzBE[1], xyzBE[2]});
216 aPoint.SetDiameter(G4VMarker::SizeType::world, pointSize);
217 aPoint.SetFillStyle(G4VMarker::FillStyle::filled);
218 sceneHandler.AddPrimitive(aPoint);
222 boxPoly.SetVisAttributes(boxZeroVisE);
223 sceneHandler.AddPrimitive(boxPoly);
225 sceneHandler.EndPrimitives();
228 G4ThreeVector startPoint = midPoint - 0.5*arrowLength*unitE;
229 G4ThreeVector endPoint = midPoint + 0.5*arrowLength*unitE;
230 G4double normalisedValue = eMag / maxFieldE;
231 if (!std::isfinite(normalisedValue))
232 {normalisedValue = 0.0;}
233 G4Colour arrowColour = eFieldColour.GetValue(normalisedValue);
236 G4Colour boxColour(arrowColour.GetRed(), arrowColour.GetGreen(), arrowColour.GetBlue(), boxAlpha);
237 sceneHandler.BeginPrimitives();
238 boxPoly.SetVisAttributes(G4VisAttributes(boxColour));
239 sceneHandler.AddPrimitive(boxPoly);
240 sceneHandler.EndPrimitives();
244 std::stringstream ss;
245 ss << std::setw(2) << std::fixed << xyzBE[0] <<
"_" << xyzBE[1] <<
"_" << xyzBE[2];
246 std::string arrowName = ss.str();
247 G4ArrowModel FArrow(startPoint.x(), startPoint.y(), startPoint.z(),
248 endPoint.x(), endPoint.y(), endPoint.z(),
249 arrowWidth, arrowColour, arrowName);
250 FArrow.DescribeYourselfTo(sceneHandler);
262 G4double result = std::min({
QIL(query->xInfo),
QIL(query->yInfo),
QIL(query->zInfo)});
269 {
return std::numeric_limits<G4double>::max();}
271 G4double step = std::abs((qi.max - qi.min) / ((G4double)qi.n));
277 G4double minSize = std::numeric_limits<double>::max();
278 G4ThreeVector result(0,0,0);
279 std::array<G4bool,3> nonZero = {
false,
false,
false};
280 std::array<const BDSFieldQueryInfo::QueryDimensionInfo,3> dims = {qi->xInfo, qi->yInfo, qi->zInfo};
283 for (G4int i = 0; i < 3; i++)
285 nonZero[i] = dims[i].n > 1;
289 result[i] = dims[i].StepSize();
290 minSize = std::min(minSize, result[i]);
296 {thinSize = 0.1 * minSize;}
298 {thinSize = 1e-2 * minSize;}
299 for (G4int i = 0; i < 3; i++)
300 {result[i] = nonZero[i] ? result[i]*0.5 : thinSize;}
Colour scale based on magma colour map.
Colour scale based on viridis colour map.
Wrapper class of BDSFieldQuery that accumulates values for visualisation.
virtual void QueryField(const BDSFieldQueryInfo *query)
Reserve the vector size we need first, then proceed as normal. Avoids big copies.
virtual void CleanUp()
Call the base class method and also clear the vector of values in this class.
Holder class for all information required for a field query.
G4double CalculateArrowLength(const BDSFieldQueryInfo *query) const
Return the minimum of 0.8 x the step length in x,y,z.
const std::vector< BDSFieldQueryInfo * > queries
Cache of queries.
G4ThreeVector BoxHalfSize(const BDSFieldQueryInfo *qi) const
G4double QIL(const BDSFieldQueryInfo::QueryDimensionInfo &qi) const
Query Info Length. For one dimension, return 0.8 x step size.
G4VisAttributes * pointVisB
The vis attributes for a zero-field marker.
G4VisAttributes * pointVisE
The vis attributes for a zero-field marker.
static G4int instanceCounter