BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSBunchUserFile.cc
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*/
19#include "BDSBunchUserFile.hh"
20#include "BDSDebug.hh"
21#include "BDSException.hh"
22#include "BDSGlobalConstants.hh"
23#include "BDSIonDefinition.hh"
24#include "BDSParticleDefinition.hh"
25#include "BDSUtilities.hh"
26#include "BDSWarning.hh"
27
28#include "parser/beam.h"
29
30#include "G4IonTable.hh"
31#include "G4ParticleDefinition.hh"
32#include "G4ParticleTable.hh"
33#include "G4String.hh"
34
35#ifdef USE_GZSTREAM
36#include "src-external/gzstream/gzstream.h"
37#endif
38
39#include <algorithm>
40#include <fstream>
41#include <regex>
42#include <set>
43#include <string>
44#include <sstream>
45#include <vector>
46
47template <class T>
49 BDSBunchFileBased("userfile"),
50 nlinesIgnore(0),
51 nlinesSkip(0),
52 nLinesValidData(0),
53 particleMass(0),
54 lineCounter(0),
55 printedOutFirstTime(false),
56 anEnergyCoordinateInUse(false),
57 changingParticleType(false),
58 endOfFileReached(false),
59 matchDistrFileLength(false)
60{
61 ffact = BDSGlobalConstants::Instance()->FFact();
62 comment = std::regex("^\\s*\\#|\\!.*");
63}
64
65template<class T>
67{
69 if (distrFile.empty())
70 {throw BDSException("BDSBunchUserFile::CheckParameters", "No input file specified for userfile distribution");}
71}
72
73template<class T>
75{
76 CloseBunchFile();
77}
78
79template<class T>
81{
82 endOfFileReached = false;
83 if (!printedOutFirstTime)
84 {
85 G4cout << "BDSBunchUserFile::OpenBunchFile> opening " << distrFilePath << G4endl;
86 printedOutFirstTime = true;
87 }
88 lineCounter = 0;
89 InputBunchFile.open(distrFilePath);
90 if (!InputBunchFile.good())
91 {throw BDSException("BDSBunchUserFile::OpenBunchFile>", "Cannot open bunch file " + distrFilePath);}
92}
93
94template<class T>
96{
97 InputBunchFile.clear(); // igzstream doesn't reset eof flags when closing - do manually
98 InputBunchFile.close();
99 endOfFileReached = true;
100}
101
102template<class T>
104{
106
107 std::regex colons(":");
108 std::vector<std::string> results;
109 std::sregex_token_iterator iter(bunchFormat.begin(), bunchFormat.end(), colons, -1);
110 std::sregex_token_iterator end;
111 for (; iter != end; ++iter)
112 {
113 std::string res = (*iter).str();
114 results.push_back(res);
115 }
116 for (auto const& token : results)
117 {
118 G4String vari = "";
119 G4String rest = "";
120 if(token.substr(0,1)=="E" || token.substr(0,1)=="P")
121 {
123 if (token.substr(0,2)=="Ek")
124 {//Kinetic energy (longer name).
125 vari = token.substr(0,2);
126 rest = token.substr(2);
127 }
128 else
129 {
130 vari = token.substr(0,1);
131 rest = token.substr(1);
132 }
133 CheckAndParseUnits(vari, rest, BDS::ParseEnergyUnit);
134 }
135 else if (token.substr(0,1)=="t")
136 {
137 vari = token.substr(0, 1);
138 rest = token.substr(1);
139 CheckAndParseUnits(vari, rest, BDS::ParseTimeUnit);
140 }
141 else if ( ( token.substr(0,1)=="x" && token.substr(1,1)!="p" ) ||
142 ( token.substr(0,1)=="y" && token.substr(1,1)!="p" ) ||
143 ( token.substr(0,1)=="z" && token.substr(1,1)!="p" ) )
144 {
145 vari = token.substr(0,1);
146 rest = token.substr(1);
147 CheckAndParseUnits(vari, rest, BDS::ParseLengthUnit);
148 }
149 else if ( (token.substr(0,2)=="xp") ||
150 (token.substr(0,2)=="yp") ||
151 (token.substr(0,2)=="zp") )
152 {
153 vari = token.substr(0,2);
154 rest = token.substr(2);
155 CheckAndParseUnits(vari, rest, BDS::ParseAngleUnit);
156 }
157 else if (token.substr(0, 1) == "S")
158 {
159 vari = token.substr(0, 1);
160 rest = token.substr(1);
161 CheckAndParseUnits(vari, rest, BDS::ParseLengthUnit);
162 useCurvilinear = true;
163 }
164 else if(token.substr(0,5)=="pdgid")
165 {
167 sd.name="pdgid";
168 sd.unit=1;
169 fields.push_back(sd);
170 }
171 else if(token.substr(0,1)=="w")
172 {
173 sd.name="weight";
174 sd.unit=1;
175 fields.push_back(sd);
176 }
177 else if (token.substr(0,1)=="-")
178 {
179 // skip
180 sd.name="skip";
181 sd.unit=1;
182 fields.push_back(sd);
183 }
184 else
185 {
186 G4String message = "Cannot determine bunch data format. Failed at token: " + token;
187 throw BDSException(__METHOD_NAME__, message);
188 }
189 }
190
191 std::set<G4String> energyKeys = {"E", "Ek", "P"};
192 CheckConflictingParameters(energyKeys);
193
194 std::set<G4String> zKeys = {"z", "S"};
196}
197
198template <typename T>
199void BDSBunchUserFile<T>::CheckConflictingParameters(const std::set<G4String>& s) const
200{
201 auto inSet = [=](const G4String& v){return s.find(v) != s.end();};
202 auto inSetD = [=](const BDSBunchUserFile::Doublet& v){return inSet(v.name);};
203 int count = std::count_if(fields.begin(), fields.end(), inSetD);
204 if (count > 1)
205 {
206 G4String message = "More than one of the following set in user file columns (\"distrFileFormat\") ";
207 for (const auto& st : s)
208 {
209 message += st;
210 message += ", ";
211 }
212 message += "\nPossibly conflicting information. Ensure only one by skipping others with \"-\" symbol";
213 throw BDSException("BDSBunchUserFile::CheckConflictingParameters>", message);
214 }
215}
216
217template <typename T>
218template <typename U>
219void BDSBunchUserFile<T>::CheckAndParseUnits(const G4String& uName,
220 const G4String& rest,
221 U unitParser)
222{
224 if (rest.empty())
225 {
226 sd.name = uName;
227 sd.unit = 1.0;
228 fields.push_back(sd);
229 return;
230 }
231
232 G4int pos1 = (G4int)rest.find("[");
233 G4int pos2 = (G4int)rest.find("]");
234 if (pos1 < 0 || pos2 < 0)
235 {
236 G4String message = "Missing bracket [] in units of \"distrFileFormat\"\n";
237 message += "variable : \"" + uName + "\" and unit \"" + rest + "\"";
238 throw BDSException("BDSBunchUserFile::CheckAndParseUnits>", message);
239 }
240 else
241 {
242 G4String fmt = rest.substr(pos1 + 1, pos2 - 1);
243 sd.name = uName;
244 sd.unit = unitParser(fmt);
245 fields.push_back(sd);
246 }
247}
248
249template<class T>
250void BDSBunchUserFile<T>::skip(std::stringstream& ss, G4int nValues)
251{
252 G4double dummyValue;
253 for (G4int i = 0; i < nValues; i++)
254 {ReadValue(ss, dummyValue);}
255}
256
257template<class T>
259{
261 {
262 if (usualPrintOut)
263 {G4cout << "BDSBunchUserFile> ignoring " << nlinesIgnore << " lines" << G4endl;}
264 std::string line;
265 for (G4int i = 0; i < (G4int)nlinesIgnore; i++)
266 {
267 std::getline(InputBunchFile, line);
268 lineCounter++;
269 // we must check explicitly if we've gone past the end of the file
270 if (InputBunchFile.eof())
271 {
272 G4String msg = "end of file reached after line " + std::to_string(lineCounter);
273 msg += " before nlinesIgnore ("+std::to_string(nlinesIgnore)+") was reached.";
274 throw BDSException("BDSBunchUserFile::SkipNLinesIgnoreIntoFile>", msg);
275 }
276 }
277 }
278}
279
280template<class T>
281void BDSBunchUserFile<T>::SkipNLinesSkip(G4bool usualPrintOut)
282{
284 {
285 if (usualPrintOut)
286 {G4cout << "BDSBunchUserFile> skipping " << nlinesSkip << " valid lines" << G4endl;}
287
288 // We can read into the file safely without checking eof() because we know from earlier
289 // counting of the number of valid lines in the file that nlinesSkip is not beyond the
290 // end of the file (including nlinesIgnore).
291 std::string line;
292 G4int nLinesValidRead = 0;
293 while (nLinesValidRead < nlinesSkip)
294 {
295 std::getline(InputBunchFile, line);
296 if (SkippableLine(line))
297 {continue;}
298 nLinesValidRead++;
299 }
300 IncrementNEventsInFileSkipped((unsigned long long int)nlinesSkip);
301 }
302}
303
304template<class T>
306 const GMAD::Beam& beam,
307 const BDSBunchType& distrType,
308 G4Transform3D beamlineTransformIn,
309 const G4double beamlineSIn)
310{
311 BDSBunchFileBased::SetOptions(beamParticle, beam, distrType, beamlineTransformIn, beamlineSIn);
312 particleMass = beamParticle->Mass();
313 distrFile = beam.distrFile;
316 nlinesIgnore = (G4long)beam.nlinesIgnore;
317 nlinesSkip = (G4long)beam.nlinesSkip;
318 matchDistrFileLength = beam.distrFileMatchLength;
320}
321
322template<class T>
323G4bool BDSBunchUserFile<T>::SkippableLine(const std::string& line) const
324{
325 return std::all_of(line.begin(), line.end(), isspace) || std::regex_search(line, comment);
326}
327
328template<class T>
330{
333
334 std::string line;
335 G4long nLinesValid = 0;
336 while ( std::getline(InputBunchFile, line) )
337 {
338 if (SkippableLine(line))
339 {continue;} // don't count it
340 nLinesValid++;
341 }
343 return nLinesValid;
344}
345
346template<class T>
348{
349 nLinesValidData = CountNLinesValidDataInFile();
350
351 if (nLinesValidData < nlinesSkip)
352 {
353 G4String msg = "nlinesSkip is greater than the number of valid lines (" + std::to_string(nLinesValidData);
354 msg += ") in the user file \"" + distrFilePath + "\"";
355 throw BDSException("BDSBunchUserFile::Initialise>", msg);
356 }
357
359 G4bool nGenerateHasBeenSet = g->NGenerateSet();
360 G4int nEventsPerLoop = (G4int)(nLinesValidData - nlinesSkip);
361 nEventsInFile = nEventsPerLoop;
362 G4int nAvailable = nEventsPerLoop * distrFileLoopNTimes;
363 G4int nGenerate = g->NGenerate();
364 if (matchDistrFileLength)
365 {
366 if (!nGenerateHasBeenSet)
367 {
368 g->SetNumberToGenerate(nAvailable);
369 G4cout << "BDSBunchUserFile::Initialise> distrFileMatchLength is true -> simulating " << nEventsPerLoop << " events";
370 if (distrFileLoopNTimes > 1)
371 {G4cout << " " << distrFileLoopNTimes << " times";}
372 G4cout << G4endl;
373 if (g->Recreate())
374 {// have to do this now before the primary generator action is called already in the run
375 // we need to start the run with a number of events so we need to calculate here
376 G4int nEventsRemaining = nAvailable - g->StartFromEvent();
377 g->SetNumberToGenerate(nEventsRemaining);
378 G4cout << "BDSBunchUserFile::Initialise> distrFileMatchLength + recreation -> simulate the "
379 << nEventsRemaining << " lines left given startFromEvent including possible looping" << G4endl;
380 }
381 }
382 else
383 {// e.g. if recreating a lower number of events - match is on; but ngenerate is lower - must obey
384 G4cout << "BDSBunchUserFile::Initialise> matchDistrFileLength has been requested "
385 << "but ngenerate has been specified -> use ngenerate" << G4endl;
386 // note we don't need to take care of a recreation offset - this is done later in primary generator action
387 if (nGenerate > nAvailable)
388 {
389 G4String msg = "ngenerate (" + std::to_string(nGenerate) + ") is greater than the number of valid lines (";
390 msg += std::to_string(nLinesValidData) + ") and distrFileMatchLength is on.\nChange ngenerate to <= # lines";
391 msg += ", or don't specifcy ngenerate.\n";
392 msg += "This includes nlinesSkip.";
393 throw BDSException("BDSBunchUserFile::Initialise>", msg);
394 }
395 }
396 }
397 else
398 {
399 if ((nGenerate > nEventsPerLoop) && !distrFileLoop)
400 {
401 G4String msg = "ngenerate (" + std::to_string(nGenerate) + ") is greater than the number of valid lines (";
402 msg += std::to_string(nLinesValidData) + ") but distrFileLoop is false in the beam command";
403 throw BDSException("BDSBunchUserFile::Initialise>", msg);
404 }
405 }
409}
410
411template<class T>
413{
414 // If the end of the file is reached go back to the beginning of the file.
415 // this re-reads the same file again - must always print warning
416 G4cout << "BDSBunchUserFile> End of file reached." << G4endl;
418 if (distrFileLoop)
419 {
420 G4cout << "BDSBunchUserFile> Returning to beginning of file (including nlinesIgnore & nlinesSkip) for next event." << G4endl;
423 SkipNLinesSkip(false);
424 }
425 else
426 {throw BDSException(__METHOD_NAME__, "distrFileLoop off but requesting another set of coordinates.");}
427}
428
429template<class T>
431{
433 G4cout << "BDSBunchUserFile::RecreateAdvanceToEvent> Advancing file to event: " << eventOffset << G4endl;
434 G4int nEventsPerLoop = nLinesValidData - nlinesSkip;
435 G4int nAvailable = nEventsPerLoop * distrFileLoopNTimes;
436 G4int nEventsRemaining = nAvailable - eventOffset;
437 if (eventOffset > nEventsPerLoop)
438 {
439 if (distrFileLoop)
440 {
441 // we're recreating a file where we looped on the data - don't repeatedly read - just read up to
442 // the right point as if we'd looped the file
443 G4int nToRoll = eventOffset % nEventsPerLoop;
444 eventOffset = nToRoll;
445 }
446 else
447 {
448 G4String msg = "eventOffset (" + std::to_string(eventOffset) + ") is greater than the number of valid data lines in this file.\n";
449 msg += "This includes nlinesSkip.";
450 throw BDSException("BDSBunchUserFile::RecreateAdvanceToEvent>", msg);
451 }
452 }
453
454 G4int nGenerate = BDSGlobalConstants::Instance()->NGenerate();
455 if (nGenerate > nEventsRemaining && !distrFileLoop)
456 {
457 G4String msg = "ngenerate (" + std::to_string(nGenerate) + ") requested in recreate mode is greater than number\n";
458 msg += "of remaining valid lines in file (" + std::to_string(nEventsRemaining) + ") and distrFileLoop is turned off.";
459 throw BDSException("BDSBunchUserFile>", msg);
460 }
461 // note we cannot update ngenerate here as we're already being called from the primary
462 // generator action in the start of the event after BeamOn(nEvents) has been called
463 // therefore this adjustment for recreation + match is done earlier in this class
464
465 // we should now be completely safe to read into the file ignoring comment lines and
466 // without checking eof()
467 std::string line;
468 for (G4int i = 0; i < eventOffset; i++)
469 {
470 std::getline(InputBunchFile, line);
471 if (SkippableLine(line))
472 {continue;}
473 lineCounter++;
474 }
475}
476
477template<class T>
479{
480 // no looping - just read one particle from file
481 return GetNextParticle();
482}
483
484template<class T>
486{
487 if (InputBunchFile.eof() || InputBunchFile.fail())
488 {EndOfFileAction();}
489
490 G4double E = 0, Ek = 0, P = 0, x = 0, y = 0, z = 0, xp = 0, yp = 0, zp = 0, t = 0;
491 G4double weight = 1;
492 G4int type = 0;
493
494 G4bool zpdef = false; //keeps record whether zp has been read from file
495 G4bool tdef = false; //keeps record whether t has been read from file
496
497 // flag whether we're going to update the particle definition
498 G4bool updateParticleDefinition = false;
499
500 // read a whole line at a time for safety - no partially read lines
501 std::string line;
502 std::getline(InputBunchFile, line);
503 lineCounter++;
504
505 // skip empty lines and comment lines (starting with # or !)
506 G4bool lineIsBad = true;
507 while (lineIsBad)
508 {
509 if (SkippableLine(line))
510 {
511 if (InputBunchFile.eof() || InputBunchFile.fail())
512 {EndOfFileAction();}
513 else
514 {
515 std::getline(InputBunchFile, line);
516 lineCounter++;
517 continue;
518 }
519 }
520 else
521 {lineIsBad = false;}
522 }
523
524 // no check the line has the right number of 'words' in it ->
525 // split line on white space - doesn't inspect words themselves
526 // checks number of words, ie number of columns is correct
527 std::vector<std::string> results;
528 std::regex wspace("\\s+"); // any whitepsace
529 // -1 here makes it point to the suffix, ie the word rather than the wspace
530 std::sregex_token_iterator iter(line.begin(), line.end(), wspace, -1);
531 std::sregex_token_iterator end;
532 for (; iter != end; ++iter)
533 {
534 std::string res = (*iter).str();
535 results.push_back(res);
536 }
537
538 if (results.size() < fields.size())
539 {// ensure enough columns
540 std::string message = "Invalid line at line " + std::to_string(lineCounter) +
541 ". Expected " + std::to_string(fields.size()) +
542 " columns , but got " + std::to_string(results.size()) +
543 ".";
544 throw BDSException(__METHOD_NAME__, message);
545 }
546
547 std::stringstream ss(line);
548 for (auto it=fields.begin();it!=fields.end();it++)
549 {
550 if(it->name=="skip")
551 {double dummy; ReadValue(ss, dummy);}
552 else if(it->name=="Ek")
553 {ReadValue(ss, Ek); Ek *= (CLHEP::GeV * it->unit);}
554 else if(it->name=="E")
555 {ReadValue(ss, E); E *= (CLHEP::GeV * it->unit);}
556 else if(it->name=="P")
557 {ReadValue(ss, P); P *= (CLHEP::GeV * it->unit);}
558 else if(it->name=="t")
559 {ReadValue(ss, t); t *= (CLHEP::s * it->unit); tdef = true;}
560 else if(it->name=="x")
561 {ReadValue(ss, x); x *= (CLHEP::m * it->unit);}
562 else if(it->name=="y")
563 {ReadValue(ss, y); y *= (CLHEP::m * it->unit);}
564 else if(it->name=="z")
565 {ReadValue(ss, z); z *= (CLHEP::m * it->unit);}
566 else if(it->name=="xp") {ReadValue(ss, xp); xp *= ( CLHEP::radian * it->unit );}
567 else if(it->name=="yp") {ReadValue(ss, yp); yp *= ( CLHEP::radian * it->unit );}
568 else if(it->name=="zp") {ReadValue(ss, zp); zp *= ( CLHEP::radian * it->unit ); zpdef = true;}
569 else if(it->name=="pdgid")
570 {// particle type
571 ReadValue(ss, type);
572 updateParticleDefinition = true; // update particle definition after finished reading line
573 }
574 else if (it->name == "S")
575 {
576 ReadValue(ss, z);
577 z *= CLHEP::m * it->unit;
578 }
579 else if(it->name=="weight")
580 {ReadValue(ss, weight);}
581 }
582
583 // coordinate checks
585 {
587 {// no checking here that only one variable is set - must ensure in this class
590 }
591 else // we must update this call of this function to ensure valid particle definition
592 {updateParticleDefinition = true;}
593 }
594 else // If energy isn't specified, use the central beam energy
595 {E = E0;}
596
597 // compute zp from xp and yp if it hasn't been read from file
598 xp += Xp0;
599 yp += Yp0;
600 if (!zpdef)
601 {zp = CalculateZp(xp,yp,1);}
602 // compute t from z if it hasn't been read from file
603 if (!tdef)
604 {t=0;} //TBC
605
606 if (updateParticleDefinition || changingParticleType)
607 {
608 // type is an int so FindParticle(int) is used here
609 G4ParticleTable* particleTable = G4ParticleTable::GetParticleTable();
610 G4ParticleDefinition* particleDef = nullptr;
611 BDSIonDefinition* ionDef = nullptr;
612 if (type < 1e9) // not a pdg ion
613 {particleDef = particleTable->FindParticle(type);}
614 else
615 {
616 G4IonTable* ionTable = particleTable->GetIonTable();
617 G4int ionA, ionZ, ionLevel;
618 G4double ionE;
619 G4IonTable::GetNucleusByEncoding(type, ionZ, ionA, ionE, ionLevel);
620 ionDef = new BDSIonDefinition(ionA, ionZ, ionZ);
621 particleDef = ionTable->GetIon(ionDef->Z(), ionDef->A(), ionDef->ExcitationEnergy());
622 }
623
624 if (!particleDef)
625 {throw BDSException("BDSBunchUserFile> Particle \"" + std::to_string(type) + "\" not found");}
626 // Wrap in our class that calculates momentum and kinetic energy.
627 // Requires that one of E, Ek, P be non-zero (only one).
628 delete particleDefinition;
629 try
630 {
631 particleDefinition = new BDSParticleDefinition(particleDef, E, Ek, P, ffact, ionDef);
634 }
635 catch (const BDSException& e)
636 {// if we throw an exception the object is invalid for the delete on the next loop
637 particleDefinition = nullptr; // reset back to nullptr for safe delete
638 throw e; // rethrow
639 }
640 }
641
642 return BDSParticleCoordsFull(X0+x,Y0+y,Z0+z,xp,yp,zp,t,z,E,weight);
643}
644
645template <class T>
646template <typename Type>
647void BDSBunchUserFile<T>::ReadValue(std::stringstream& stream, Type& value)
648{
649 stream >> value;
650}
651
653
654#ifdef USE_GZSTREAM
655template class BDSBunchUserFile<igzstream>;
656#endif
657
658namespace BDS
659{
660G4double ParseEnergyUnit(const G4String& fmt)
661{
662 G4double unit= 1.0;
663 if (fmt == "TeV")
664 {unit = 1e3;}
665 else if (fmt == "GeV")
666 {unit = 1.0 ;}
667 else if (fmt == "MeV")
668 {unit = 1e-3;}
669 else if (fmt == "KeV")
670 {unit = 1e-6;}
671 else if (fmt == "eV")
672 {unit = 1e-9;}
673 else
674 {throw BDSException(__METHOD_NAME__, "Unrecognised energy unit! " +fmt);}
675 return unit;
676}
677
678G4double ParseLengthUnit(const G4String& fmt)
679{
680 G4double unit = 1.0;
681 if (fmt == "m")
682 {unit = 1;}
683 else if (fmt == "cm")
684 {unit = 1e-2;}
685 else if (fmt == "mm")
686 {unit = 1e-3;}
687 else if (fmt == "mum" || fmt == "um")
688 {unit = 1e-6;}
689 else if (fmt == "nm")
690 {unit = 1e-9;}
691 else
692 {throw BDSException(__METHOD_NAME__, "Unrecognised length unit! " + fmt);}
693 return unit;
694}
695
696G4double ParseAngleUnit(const G4String& fmt)
697{
698 G4double unit = 1.0;
699 if (fmt == "rad")
700 {unit = 1;}
701 else if (fmt == "mrad")
702 {unit = 1e-3;}
703 else if (fmt == "murad" || fmt == "urad")
704 {unit = 1e-6;}
705 else if (fmt == "nrad")
706 {unit=1e-9;}
707 else
708 {throw BDSException(__METHOD_NAME__, "Unrecognised angle unit! " + fmt);}
709 return unit;
710}
711
712G4double ParseTimeUnit(const G4String& fmt)
713{
714 G4double unit = 1.0;
715 if (fmt == "s")
716 {unit = 1;}
717 else if (fmt == "ms")
718 {unit = 1.e-3;}
719 else if (fmt == "mus" || fmt == "us")
720 {unit = 1e-6;}
721 else if (fmt == "ns")
722 {unit = 1e-9;}
723 else if (fmt == "mm/c")
724 {unit = (CLHEP::mm/CLHEP::c_light)/CLHEP::s;}
725 else if (fmt == "nm/c")
726 {unit = (CLHEP::nm/CLHEP::c_light)/CLHEP::s;}
727 else
728 {throw BDSException(__METHOD_NAME__, "Unrecognised time unit! " + fmt);}
729 return unit;
730}
731
732}
An intermediate layer for any bunch classes that are file based.
unsigned long long int nEventsInFile
The number of entries in the file loaded.
virtual void SetOptions(const BDSParticleDefinition *beamParticle, const GMAD::Beam &beam, const BDSBunchType &distrType, G4Transform3D beamlineTransformIn=G4Transform3D::Identity, const G4double beamlineS=0)
Pull out the relevant options and then pass through to BDSBunch::SetOptions().
A bunch distribution that reads a user specified column file.
void CheckConflictingParameters(const std::set< G4String > &s) const
Check conflicting columns aren't specified in file, e.g. P and Ek. Throw exception if wrong.
virtual void RecreateAdvanceToEvent(G4int eventOffset)
G4double ffact
Cache of flip factor from global constants.
void ReadValue(std::stringstream &stream, Type &value)
G4int lineCounter
Line counter.
G4long CountNLinesValidDataInFile()
void SkipNLinesSkip(G4bool usualPrintOut=true)
T InputBunchFile
The file handler. Templated as could be std::ifstream or igzstream for example.
G4bool SkippableLine(const std::string &line) const
Return true if a line is all whitespace or is commented out (starts with '#').
virtual BDSParticleCoordsFullGlobal GetNextParticleValid(G4int maxTries)
G4bool changingParticleType
Whether the particle type is a column.
void skip(std::stringstream &stream, G4int nvalues)
void SkipNLinesIgnoreIntoFile(G4bool usualPrintOut=true)
void ParseFileFormat()
Parse the column tokens and units factors.
G4long nlinesIgnore
Number of lines that will be ignored at the start the file.
virtual void SetOptions(const BDSParticleDefinition *beamParticle, const GMAD::Beam &beam, const BDSBunchType &distrType, G4Transform3D beamlineTransformIn=G4Transform3D::Identity, const G4double beamlineS=0)
Pull out the relevant options and then pass through to BDSBunch::SetOptions().
G4String distrFilePath
Bunch file including absolute path.
G4String bunchFormat
Format of the file.
G4long nlinesSkip
Number of lines that will be skipped after the nlinesIgnore.
G4double particleMass
Cache of nominal beam particle mass.
std::list< Doublet > fields
List of variables to parse on each line.
void CloseBunchFile()
Close the file handler.
virtual void CheckParameters()
virtual BDSParticleCoordsFull GetNextParticleLocal()
Get the next particle.
virtual void Initialise()
Open the file and skip lines.
G4String distrFile
Bunch file.
void OpenBunchFile()
Open the file and check it's open.
G4bool anEnergyCoordinateInUse
Whether Et, Ek or P are in the columns.
G4double Yp0
Centre of distributions.
Definition: BDSBunch.hh:177
virtual void RecreateAdvanceToEvent(G4int eventOffset)
Definition: BDSBunch.cc:280
G4bool particleDefinitionHasBeenUpdated
Definition: BDSBunch.hh:203
G4bool useCurvilinear
Whether to ignore z and use s and transform for curvilinear coordinates.
Definition: BDSBunch.hh:196
G4double Z0
Centre of distributions.
Definition: BDSBunch.hh:173
G4double X0
Centre of distributions.
Definition: BDSBunch.hh:171
G4double Xp0
Centre of distributions.
Definition: BDSBunch.hh:176
G4double E0
Centre of distributions.
Definition: BDSBunch.hh:179
G4double Y0
Centre of distributions.
Definition: BDSBunch.hh:172
BDSParticleCoordsFullGlobal GetNextParticle()
Definition: BDSBunch.cc:262
static G4double CalculateZp(G4double xp, G4double yp, G4double Zp0)
Calculate zp safely based on other components.
Definition: BDSBunch.cc:367
BDSParticleDefinition * particleDefinition
Particle definition for bunch - this class owns it.
Definition: BDSBunch.hh:199
virtual void CheckParameters()
Definition: BDSBunch.cc:223
General exception with possible name of object and message.
Definition: BDSException.hh:35
static BDSGlobalConstants * Instance()
Access method.
Class to parse an ion particle definition.
G4int A() const
Accessor.
G4double ExcitationEnergy() const
Accessor.
G4int Z() const
Accessor.
A set of particle coordinates in both local and global.
A set of particle coordinates including energy and weight.
Wrapper for particle definition.
G4double Mass() const
Accessor.
void SetEnergies(G4double totalEnergyIn, G4double kineticEnergyIn, G4double momentumIn)
G4double TotalEnergy() const
Accessor.
Improve type-safety of native enum data type in C++.
std::string distrFileFormat
beam parameters
Definition: beamBase.h:53
std::string distrFile
beam parameters
Definition: beamBase.h:52
int nlinesIgnore
Ignore first lines in the input bunch file.
Definition: beamBase.h:60
int nlinesSkip
Number of event lines to skip after the ignore lines.
Definition: beamBase.h:61
bool distrFileMatchLength
beam parameters
Definition: beamBase.h:55
Beam class.
Definition: beam.h:44
Return either G4Tubs or G4CutTubs depending on flat face.
G4String GetFullPath(G4String filename, bool excludeNameFromPath=false, bool useCWDForPrefix=false)
G4double ParseTimeUnit(const G4String &fmt)
G4double ParseEnergyUnit(const G4String &fmt)
G4double ParseAngleUnit(const G4String &fmt)
G4bool IsFinite(G4double value, G4double tolerance=std::numeric_limits< double >::epsilon())
G4double ParseLengthUnit(const G4String &fmt)
Struct for name and unit pair.
G4double unit
relative to SI units, i.e. mm=0.001 etc.