EGSnrc C++ class library  Report PIRS-898 (2021)
Iwan Kawrakow, Ernesto Mainegra-Hing, Frederic Tessier, Reid Townson and Blake Walters
egs_dynamic_source.cpp
Go to the documentation of this file.
1 /*
2 ###############################################################################
3 #
4 # EGSnrc egs++ dynamic source
5 # Copyright (C) 2015 National Research Council Canada
6 #
7 # This file is part of EGSnrc.
8 #
9 # EGSnrc is free software: you can redistribute it and/or modify it under
10 # the terms of the GNU Affero General Public License as published by the
11 # Free Software Foundation, either version 3 of the License, or (at your
12 # option) any later version.
13 #
14 # EGSnrc is distributed in the hope that it will be useful, but WITHOUT ANY
15 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
16 # FOR A PARTICULAR PURPOSE. See the GNU Affero General Public License for
17 # more details.
18 #
19 # You should have received a copy of the GNU Affero General Public License
20 # along with EGSnrc. If not, see <http://www.gnu.org/licenses/>.
21 #
22 ###############################################################################
23 #
24 # Author: Blake Walters, 2017
25 #
26 # Contributors: Alexandre Demelo
27 #
28 ###############################################################################
29 */
30 
31 
37 #include "egs_dynamic_source.h"
38 #include "egs_input.h"
39 
41  EGS_ObjectFactory *f) : EGS_BaseSource(input,f), source(0), valid(true) {
42  EGS_Input *isource = input->takeInputItem("source",false);
43  if (isource) {
45  delete isource;
46  }
47  if (!source) {
48  string sname;
49  int err = input->getInput("source name",sname);
50  if (err)
51  egsWarning("EGS_DynamicSource: missing/wrong inline source "
52  "definition and missing wrong 'source name' input\n");
53  else {
55  if (!source) egsWarning("EGS_DynamicSource: a source named %s"
56  " does not exist\n",sname.c_str());
57  }
58  }
59  //now read inputs relevant to dynamic source
60  //see if user wants to synchronize source with time read from
61  //iaea phsp or beam simulation source
62  vector<string> sync_options;
63  sync_options.push_back("no");
64  sync_options.push_back("yes");
65  sync = input->getInput("synchronize motion",sync_options,0);
66  if (sync && source->getObjectType()!="IAEA_PhspSource" &&
67  source->getObjectType()!="EGS_BeamSource") {
68  egsWarning("EGS_DynamicSource: source motion can only be synchronized with a source of type iaea_phsp_source or egs_beam_source.\n Will not synchronize.\n");
69  sync = false;
70  }
71 
72  //get control points
73  EGS_Input *dyninp = input->takeInputItem("motion");
74  if (dyninp) {
75  //get control points
76  ncpts=0;
77  vector<EGS_Float> point;
78  EGS_ControlPoint cpt;
79  stringstream itos;
80  int err;
81  int icpts=1;
82  itos << icpts;
83 
84  string inputTag = "control point";
85  string inputTag_backCompat = "control point " + itos.str();
86  EGS_Input *currentInput;
87 
88  while (true) {
89  currentInput = dyninp->takeInputItem(inputTag);
90  if (!currentInput || currentInput->getInput(inputTag, point)) {
91  currentInput = dyninp->takeInputItem(inputTag_backCompat);
92  if (!currentInput || currentInput->getInput(inputTag_backCompat, point)) {
93  delete currentInput;
94  break;
95  }
96  }
97  delete currentInput;
98 
99  if (point.size()!=8) {
100  egsWarning("EGS_DynamicSource: control point %i does not specify 8 values.\n",icpts);
101  valid = false;
102  }
103  else {
104  if (ncpts>0 && point[7] < cpts[ncpts-1].time) {
105  egsWarning("EGS_DynamicSource: time index of control point %i < time index of control point %i\n",icpts,ncpts);
106  valid = false;
107  }
108  else if (point[7] < 0.) {
109  egsWarning("EGS_DynamicSource: time index of control point %i < 0.0\n",icpts);
110  valid = false;
111  }
112  else {
113  ncpts++;
114  if (ncpts ==1 && point[7] > 0.0) {
115  egsWarning("EGS_DynamicSource: time index of control point 1 > 0.0. This will generate many warning messages.\n");
116  }
117  //set cpt values
118  cpt.iso = EGS_Vector(point[0],point[1],point[2]);
119  cpt.dsource = point[3];
120  cpt.theta = point[4];
121  cpt.phi = point[5];
122  cpt.phicol = point[6];
123  cpt.time = point[7];
124  //add it to the vector of control points
125  cpts.push_back(cpt);
126  icpts++;
127  itos.str("");
128  itos << icpts;
129  inputTag_backCompat = "control point " + itos.str();
130  }
131  }
132  }
133  if (ncpts<=1) {
134  egsWarning("EGS_DynamicSource: not enough or missing control points.\n");
135  valid = false;
136  }
137  if (cpts[ncpts-1].time == 0.0) {
138  egsWarning("EGS_DynamicSource: time index of last control point = 0. Something's wrong.\n");
139  valid = false;
140  }
141  else {
142  //normalize time index to max. value
143  for (int i=0; i<=ncpts-1; i++) {
144  cpts[i].time /= cpts[ncpts-1].time;
145  }
146  }
147  }
148  else {
149  egsWarning("EGS_DynamicSource: no control points input.\n");
150  valid = false;
151  }
152  setUp();
153 }
154 
155 void EGS_DynamicSource::setUp() {
156  //most setup done in constructor
157  otype="EGS_DynamicSource";
158  if (!isValid()) {
159  description = "Invalid dynamic source";
160  }
161  else {
162  description = "Dynamic source based on\n";
164  if (sync) {
165  description += "\n Source will be synched with time values read in (if available).";
166  }
167  }
168 }
169 
170 //actually select the rotation coordinates for the incident particle
171 int EGS_DynamicSource::getCoord(EGS_Float rand, EGS_ControlPoint &ipt) {
172  int iindex=0;
173  int i;
174  for (i=0; i<ncpts; i++) {
175  if (rand < cpts[i].time-epsilon) {
176  iindex =i;
177  break;
178  }
179  }
180  if (i==ncpts) {
181  egsWarning("EGS_DynamicSource: could not locate control point.\n");
182  return 1;
183  }
184  EGS_Float factor = (rand-cpts[iindex-1].time)/(cpts[iindex].time-cpts[iindex-1].time);
185  ipt.iso.x=cpts[iindex-1].iso.x+ (cpts[iindex].iso.x-cpts[iindex-1].iso.x)*factor;
186  ipt.iso.y=cpts[iindex-1].iso.y+ (cpts[iindex].iso.y-cpts[iindex-1].iso.y)*factor;
187  ipt.iso.z=cpts[iindex-1].iso.z+ (cpts[iindex].iso.z-cpts[iindex-1].iso.z)*factor;
188  ipt.dsource=cpts[iindex-1].dsource+ (cpts[iindex].dsource-cpts[iindex-1].dsource)*factor;
189  ipt.theta=cpts[iindex-1].theta+ (cpts[iindex].theta-cpts[iindex-1].theta)*factor;
190  ipt.phi=cpts[iindex-1].phi+ (cpts[iindex].phi-cpts[iindex-1].phi)*factor;
191  ipt.phicol=cpts[iindex-1].phicol+ (cpts[iindex].phicol-cpts[iindex-1].phicol)*factor;
192  return 0;
193 };
194 
200 void EGS_DynamicSource::containsDynamic(bool &hasdynamic) {
201  hasdynamic = true;
202 }
203 
204 extern "C" {
205 
206  EGS_DYNAMIC_SOURCE_EXPORT EGS_BaseSource *createSource(EGS_Input *input,
207  EGS_ObjectFactory *f) {
208  return
209  createSourceTemplate<EGS_DynamicSource>(input,f,"dynamic source");
210  }
211 
212 }
Base source class. All particle sources must be derived from this class.
const char * getSourceDescription() const
Get a short description of this source.
static EGS_BaseSource * getSource(const string &Name)
Get a pointer to the source named Name.
string description
A short source description.
static EGS_BaseSource * createSource(EGS_Input *)
Create sources from the information pointed to by input.
EGS_BaseSource * source
The source being rotated.
EGS_DynamicSource(EGS_BaseSource *Source, vector< EGS_ControlPoint > cpts, const string &Name="", EGS_ObjectFactory *f=0)
Construct a dynamic source using Source as the source and cpts as the control points....
void containsDynamic(bool &hasdynamic)
Check if the simulation source contains time indices.
A class for storing information in a tree-like structure of key-value pairs. This class is used throu...
Definition: egs_input.h:182
EGS_Input * takeInputItem(const string &key, bool self=true)
Get the property named key.
Definition: egs_input.cpp:226
int getInput(const string &key, vector< string > &values) const
Assign values to an array of strings from an input identified by key.
Definition: egs_input.cpp:338
An object factory.
string otype
The object type.
const string & getObjectType() const
Get the object type.
A class representing 3D vectors.
Definition: egs_vector.h:57
A source with simulated time-varying rotations/translations.
EGS_Input class header file.
const EGS_Float epsilon
The epsilon constant for floating point comparisons.
Definition: egs_functions.h:62
EGS_InfoFunction EGS_EXPORT egsWarning
Always use this function for reporting warnings.