EGSnrc C++ class library  Report PIRS-898 (2021)
Iwan Kawrakow, Ernesto Mainegra-Hing, Frederic Tessier, Reid Townson and Blake Walters
egs_smart_envelope.h
Go to the documentation of this file.
1 /*
2 ###############################################################################
3 #
4 # EGSnrc egs++ smart envelope geometry headers
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: Iwan Kawrakow, 2008
25 #
26 # Contributors: Frederic Tessier
27 # Ernesto Mainegra-Hing
28 #
29 ###############################################################################
30 #
31 # A "smart" envelope geometry. The smartness is due to the fact that there
32 # can be only zero or one inscribed geometry per base geometry region, which
33 # makes many of the checks faster.
34 #
35 # In addition, unlike the regular envelope geometry where all inscribed
36 # geometries always must completely fit inside the base geometry, here
37 # geometries can be inscribed using logic 0 or 1.
38 #
39 # Logic 0 is as before, i.e., geometry must completely fit into the region
40 # where it is being inscribed.
41 #
42 # Logic 1 means that the inscribed geometry extends beyond the region and
43 # therefore the smart envelope also checks the base geometry in this case.
44 # This is very similar to a CD geometry except that now the space outside
45 # the inscribed geometry is still part of the geometry.
46 #
47 # Warning: not completely tested, so don't use for production runs yet.
48 #
49 ###############################################################################
50 */
51 
52 
58 #ifndef EGS_SMART_ENVELOPE_
59 #define EGS_SMART_ENVELOPE_
60 
61 #ifdef WIN32
62 
63  #ifdef BUILD_SMART_ENVELOPE_DLL
64  #define EGS_SMART_ENVELOPE_EXPORT __declspec(dllexport)
65  #else
66  #define EGS_SMART_ENVELOPE_EXPORT __declspec(dllimport)
67  #endif
68  #define EGS_SMART_ENVELOPE_LOCAL
69 
70 #else
71 
72  #ifdef HAVE_VISIBILITY
73  #define EGS_SMART_ENVELOPE_EXPORT __attribute__ ((visibility ("default")))
74  #define EGS_SMART_ENVELOPE_LOCAL __attribute__ ((visibility ("hidden")))
75  #else
76  #define EGS_SMART_ENVELOPE_EXPORT
77  #define EGS_SMART_ENVELOPE_LOCAL
78  #endif
79 
80 #endif
81 
82 
83 #include "egs_base_geometry.h"
84 #include "egs_functions.h"
85 
86 #include<vector>
87 using std::vector;
88 
134 //#include "egs_functions.h"
135 
136 struct SmartEnvelopeAux;
137 
138 class EGS_SMART_ENVELOPE_EXPORT EGS_SmartEnvelope : public EGS_BaseGeometry {
139 
140 public:
141 
143  const vector<SmartEnvelopeAux *> &fgeoms, const string &Name = "");
144 
146 
147  bool isInside(const EGS_Vector &x) {
148  return g->isInside(x);
149  };
150 
151  int isWhere(const EGS_Vector &x) {
152  int ibase = g->isWhere(x);
153  if (ibase < 0) {
154  return ibase;
155  }
156  int j = gindex[ibase];
157  if (j >= 0) {
158  int i = geometries[j]->isWhere(x);
159  if (i >= 0) {
160  ibase = local_start[j] + i;
161  }
162  }
163  return ibase;
164  };
165 
166  int inside(const EGS_Vector &x) {
167  return isWhere(x);
168  };
169 
170  bool isRealRegion(int ireg) const {
171  if (ireg < nbase) {
172  return g->isRealRegion(ireg);
173  }
174  int i = ireg - nbase;
175  int j = reg_to_inscr[i];
176  return geometries[j]->isRealRegion(ireg-local_start[j]);
177  };
178 
179  int medium(int ireg) const {
180  if (ireg < nbase) {
181  return g->medium(ireg);
182  }
183  int i = ireg - nbase;
184  int j = reg_to_inscr[i];
185  return geometries[j]->medium(ireg-local_start[j]);
186  };
187 
188  EGS_Float howfarToOutside(int ireg, const EGS_Vector &x,
189  const EGS_Vector &u) {
190  if (ireg < 0) {
191  return 0;
192  }
193  int ibase = ireg < nbase ? ireg : reg_to_base[ireg-nbase];
194  return g->howfarToOutside(ibase,x,u);
195  };
196 
197  int howfar(int ireg, const EGS_Vector &x, const EGS_Vector &u,
198  EGS_Float &t, int *newmed = 0, EGS_Vector *normal = 0) {
199  if (ireg >= 0) {
200  // inside.
201  if (ireg < nbase) {
202  // in one of the regions of the base geometry
203  // check if we hit a boundary in the base geometry.
204  // if we do, newmed and normal get set accordingly.
205  int ibase = g->howfar(ireg,x,u,t,newmed,normal);
206  int j = gindex[ireg];
207  if (j >= 0) {
208  // there is an inscribed geometry in this region
209  // => check if we will enter
210  int inscr = geometries[j]->howfar(-1,x,u,t,newmed,normal);
211  if (inscr >= 0)
212  // yes, we do.
213  {
214  return local_start[j] + inscr;
215  }
216  }
217  if (ibase == ireg || ibase < 0) {
218  return ibase;
219  }
220  // if here, we enter a new base geometry region.
221  j = gindex[ibase];
222  if (j >= 0 && itype[j]) {
223  // geometry inscribed in the new base region
224  // using inscription type 1 => check if already
225  // inside the inscribed
226  int inscr = geometries[j]->isWhere(x+u*t);
227  if (inscr >= 0) {
228  if (newmed) {
229  *newmed = geometries[j]->medium(inscr);
230  }
231  return local_start[j] + inscr;
232  }
233  }
234  return ibase;
235  }
236  // if here, we are in an inscribed geometry.
237  int i = ireg-nbase;
238  int ibase = reg_to_base[i];
239  int j = reg_to_inscr[i];
240  int ilocal = ireg-local_start[j];
241  int inew = geometries[j]->howfar(ilocal,x,u,t,newmed,normal);
242  if (itype[j]) {
243  int ibase_new = g->howfar(ibase,x,u,t,newmed,normal);
244  if (ibase_new != ibase) {
245  if (ibase_new < 0) {
246  return ibase_new;
247  }
248  j = gindex[ibase_new];
249  if (j >= 0 && itype[j]) {
250  int inscr = geometries[j]->isWhere(x+u*t);
251  if (inscr >= 0) {
252  if (newmed) {
253  *newmed = geometries[j]->medium(inscr);
254  }
255  return local_start[j] + inscr;
256  }
257  }
258  return ibase_new;
259  }
260  }
261  if (inew < 0) {
262  if (newmed) {
263  *newmed = g->medium(ibase);
264  }
265  return ibase;
266  }
267  return local_start[j] + inew;
268  }
269  // if here, we are outside the base geometry.
270  // check to see if we will enter.
271  int ibase = g->howfar(ireg,x,u,t,newmed,normal);
272  if (ibase >= 0) {
273  int j = gindex[ibase];
274  if (j >= 0 && itype[j]) {
275  int inscr = geometries[j]->isWhere(x+u*t);
276  if (inscr >= 0) {
277  if (newmed) {
278  *newmed = geometries[j]->medium(inscr);
279  }
280  return local_start[j] + inscr;
281  }
282  }
283  }
284  return ibase;
285  };
286 
287  EGS_Float hownear(int ireg, const EGS_Vector &x) {
288  if (ireg >= 0) {
289  EGS_Float tmin;
290  if (ireg < nbase) { // in one of the regions of the base geom.
291  tmin = g->hownear(ireg,x);
292  if (tmin <= 0) {
293  return tmin;
294  }
295  int j = gindex[ireg];
296  if (j >= 0) {
297  EGS_Float ti = geometries[j]->hownear(-1,x);
298  if (ti < tmin) {
299  tmin = ti;
300  }
301  }
302  return tmin;
303  }
304  int i = ireg-nbase;
305  int j = reg_to_inscr[i];
306  int ilocal = ireg-local_start[j];
307  tmin = geometries[j]->hownear(ilocal,x);
308  if (itype[j]) {
309  int ibase = reg_to_base[i];
310  EGS_Float tbase = g->hownear(ibase,x);
311  if (tbase < tmin) {
312  tmin = tbase;
313  }
314  }
315  return tmin;
316  }
317  return g->hownear(ireg,x);
318  };
319 
320  int getMaxStep() const {
321  int nstep = g->getMaxStep() + n_in;
322  for (int j=0; j<n_in; ++j) {
323  nstep += geometries[j]->getMaxStep();
324  }
325  return nstep;
326  };
327 
328  bool hasBooleanProperty(int ireg, EGS_BPType prop) const {
329  if (ireg >= 0 && ireg < nreg) {
330  if (ireg < nbase) {
331  return g->hasBooleanProperty(ireg,prop);
332  }
333  int i = ireg-nbase, j = reg_to_inscr[i];
334  return geometries[j]->hasBooleanProperty(ireg-local_start[j],prop);
335  }
336  return false;
337  };
338  void setBooleanProperty(EGS_BPType) {
339  setPropertyError("setBooleanProperty()");
340  };
341  void addBooleanProperty(int) {
342  setPropertyError("addBooleanProperty()");
343  };
344  void setBooleanProperty(EGS_BPType,int,int,int step=1) {
345  setPropertyError("setBooleanProperty()");
346  };
347  void addBooleanProperty(int,int,int,int step=1) {
348  setPropertyError("addBooleanProperty()");
349  };
350 
351  const string &getType() const {
352  return type;
353  };
354 
355  void printInfo() const;
356 
357  void setRelativeRho(int start, int end, EGS_Float rho);
358  void setRelativeRho(EGS_Input *);
359  EGS_Float getRelativeRho(int ireg) const {
360  if (ireg < 0 || ireg >= nreg) {
361  return 1;
362  }
363  if (ireg < nbase) {
364  return g->getRelativeRho(ireg);
365  }
366  int i = ireg-nbase;
367  int j = reg_to_inscr[i];
368  return geometries[j]->getRelativeRho(ireg-local_start[j]);
369  };
370 
371  void setBScaling(int start, int end, EGS_Float bf);
372  void setBScaling(EGS_Input *);
373  EGS_Float getBScaling(int ireg) const {
374  if (ireg < 0 || ireg >= nreg) {
375  return 1;
376  }
377  if (ireg < nbase) {
378  return g->getBScaling(ireg);
379  }
380  int i = ireg-nbase;
381  int j = reg_to_inscr[i];
382  return geometries[j]->getBScaling(ireg-local_start[j]);
383  };
384 
385  virtual void getLabelRegions(const string &str, vector<int> &regs);
386 
387 protected:
388 
391  int *gindex;
393  int *reg_to_base;
394  int *local_start;
395  char *itype;
396  int n_in;
397  int nbase;
398 
399  static string type;
400 
407  void setMedia(EGS_Input *,int,const int *);
408 
409 private:
410 
411  void setPropertyError(const char *funcname) {
412  egsFatal("EGS_SmartEnvelope::%s: don't use this method\n Define "
413  "properties in the constituent geometries instead\n",
414  funcname);
415  };
416 
417 
418 };
419 
420 #endif
virtual const string & getType() const =0
Get the geometry type.
virtual bool isInside(const EGS_Vector &x)=0
Is the position x inside the geometry?
virtual bool hasBooleanProperty(int ireg, EGS_BPType prop) const
Is the boolean property prop set for region ireg ?
virtual int medium(int ireg) const
Returns the medium index in region ireg.
virtual void printInfo() const
Print information about this geometry.
int * gindex
Index of inscribed geometries.
virtual EGS_Float howfarToOutside(int ireg, const EGS_Vector &x, const EGS_Vector &u)
virtual int getMaxStep() const
Returns the maximum number of steps through the geometry.
A class representing 3D vectors.
Definition: egs_vector.h:56
int * reg_to_base
Region to base region conversion.
Global egspp functions header file.
virtual int isWhere(const EGS_Vector &x)=0
In which region is poisition x?
virtual int howfar(int ireg, const EGS_Vector &x, const EGS_Vector &u, EGS_Float &t, int *newmed=0, EGS_Vector *normal=0)=0
Calculate the distance to a boundary from x along the direction u.
Base geometry class. Every geometry class must be derived from EGS_BaseGeometry.
virtual void setRelativeRho(int start, int end, EGS_Float rho)
Set the relative mass density in regions.
virtual EGS_Float getBScaling(int ireg) const
Get the B field scaling factor in region ireg.
virtual void setBScaling(int start, int end, EGS_Float bf)
Set the B field scaling factor in regions.
EGS_InfoFunction EGS_EXPORT egsFatal
Always use this function for reporting fatal errors.
virtual void getLabelRegions(const string &str, vector< int > &regs)
Get the list of all regions labeled with str.
void setMedia(EGS_Input *inp)
Set the media in the geometry from the input pointed to by inp.
virtual int inside(const EGS_Vector &x)=0
Returns the region index, if inside, or -1 if outside (obsolete)
int n_in
Number of inscribed geometries.
virtual EGS_Float hownear(int ireg, const EGS_Vector &x)=0
Calculate the distance to a boundary for position x in any direction.
int * reg_to_inscr
Region to inscribed geometry conversion.
virtual EGS_Float getRelativeRho(int ireg) const
Get the relative mass density in region ireg.
EGS_BaseGeometry * g
The envelope geometry.
int * local_start
First region for each inscribed geometry.
EGS_BaseGeometry ** geometries
The inscribed geometries.
static string type
Geometry type.
virtual bool isRealRegion(int ireg) const
Returnes true if ireg is a real region, false otherwise.
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_BaseGeometry class header file.
int nbase
Number of regions in the base geometry.
virtual void addBooleanProperty(int bit)
Add a boolean property for the entire geometry by setting the bit&#39;th bit.
virtual void setBooleanProperty(EGS_BPType prop)
Set the boolean properties of the entire geometry to prop.