Add dynamic motion for shapes using egs_dynamic_shape

HEN_HOUSE/egs++/Makefile

@@ -27,6 +27,7 @@
 #                   Frederic Tessier
 #                   Manuel Stoeckl
 #                   Max Orok
+#                   Alexandre Demelo
 #
 ###############################################################################
 
@@ -53,7 +54,8 @@ geometry_libs = egs_planes egs_cd_geometry egs_gtransformed egs_nd_geometry \
                egs_cones egs_gstack egs_prism egs_union egs_pyramid egs_conez\
                egs_space egs_elliptic_cylinders egs_smart_envelope \
                egs_vhp_geometry egs_octree egs_roundrect_cylinders \
-               egs_rz egs_autoenvelope egs_lattice egs_glib egs_mesh
+               egs_rz egs_autoenvelope egs_lattice egs_glib egs_mesh \
+               egs_dynamic_geometry
 
 source_libs = egs_collimated_source egs_isotropic_source egs_parallel_beam \
               egs_point_source egs_source_collection egs_transformed_source \
@@ -64,7 +66,7 @@ source_libs = egs_collimated_source egs_isotropic_source egs_parallel_beam \
 shape_libs = egs_circle egs_ellipse egs_extended_shape egs_gaussian_shape \
              egs_line_shape egs_polygon_shape egs_rectangle egs_shape_collection \
              egs_voxelized_shape egs_spherical_shell egs_conical_shell \
-             egs_circle_perpendicular
+             egs_circle_perpendicular egs_dynamic_shape
 
 aobject_libs = egs_track_scoring egs_dose_scoring egs_radiative_splitting \
 	       egs_phsp_scoring egs_fluence_scoring

HEN_HOUSE/egs++/ausgab_objects/egs_phsp_scoring/egs_phsp_scoring.cpp

@@ -242,7 +242,7 @@ void EGS_PhspScoring::reportResults() {
 void EGS_PhspScoring::storeParticle(EGS_I64 ncase) {
 
     //if user requested mu scoring, check if mu is available
-    if (score_mu && app->getMU() < 0) {
+    if (score_mu && app->getTimeIndex() < 0) {
         egsWarning("\nEGS_PhspScoring: User requested mu scoring, but mu is inavailable with this source.\n");
         egsWarning("Turning off mu scoring.\n");
         score_mu=false;
@@ -278,7 +278,7 @@ void EGS_PhspScoring::storeParticle(EGS_I64 ncase) {
     p_stack[phsp_index].w = app->top_p.u.z;
     p_stack[phsp_index].q = app->top_p.q;
     if (score_mu) {
-        p_stack[phsp_index].mu = app->getMU();
+        p_stack[phsp_index].mu = app->getTimeIndex();
     }
     p_stack[phsp_index++].latch = app->top_p.latch;
 

HEN_HOUSE/egs++/ausgab_objects/egs_track_scoring/egs_track_scoring.cpp

@@ -24,6 +24,7 @@
 #  Author:          Iwan Kawrakow, 2009
 #
 #  Contributors:    Georgi Gerganov
+#                   Alexandre Demelo
 #
 ###############################################################################
 */
@@ -41,7 +42,7 @@
 EGS_TrackScoring::EGS_TrackScoring(const string &Name, EGS_ObjectFactory *f) :
     EGS_AusgabObject(Name,f), m_pts(0), m_start(0), m_stop(1024), m_lastCase(-1),
     m_nScore(0), m_bufSize(16), m_score(false), m_didScore(false),
-    m_score_photons(true), m_score_electrons(true), m_score_positrons(true), m_fnExtra("") {
+    m_score_photons(true), m_score_electrons(true), m_score_positrons(true), m_fnExtra(""), m_include_time(false) {
     otype = "EGS_TrackScoring";
 }
 
@@ -77,8 +78,17 @@ void EGS_TrackScoring::setApplication(EGS_Application *App) {
         sprintf(buf,"_w%d",i_parallel);
         fname += buf;
     }
-    fname += ".ptracks";
-    m_pts = new EGS_ParticleTrackContainer(fname.c_str(),m_bufSize);
+    // if incltime is false use .ptracks. If incltime is true use the new file
+    // format .syncptracks
+    if (m_include_time) {
+        fname += ".syncptracks";
+    }
+    else {
+        fname += ".ptracks";
+    }
+    // create new particleTrackContainer using the m_include_time boolean which
+    // controls time index writting and filetype
+    m_pts = new EGS_ParticleTrackContainer(fname.c_str(),m_bufSize,m_include_time);
 
     description = "\nParticle Track Scoring (";
     description += name;
@@ -90,6 +100,8 @@ void EGS_TrackScoring::setApplication(EGS_Application *App) {
     description += m_score_electrons ? "YES\n" : "NO\n";
     description += " - Scoring positron tracks     = ";
     description += m_score_positrons ? "YES\n" : "NO\n";
+    description += " - Include time index          = ";
+    description += m_include_time ? "YES\n" : "NO\n";
     description += " - First event to score        = ";
     char buf[32];
     sprintf(buf,"%lld\n",m_start);
@@ -114,7 +126,6 @@ void EGS_TrackScoring::reportResults() {
     }
 }
 
-
 extern "C" {
 
     EGS_TRACK_SCORING_EXPORT EGS_AusgabObject *createAusgabObject(EGS_Input *input,
@@ -130,6 +141,11 @@ extern "C" {
         bool scph = input->getInput("score photons",sc_options,true);
         bool scel = input->getInput("score electrons",sc_options,true);
         bool scpo = input->getInput("score positrons",sc_options,true);
+
+        // include time index lets the program know whether to write the time
+        // index to the tracks file and determines the filetype (ptracks or
+        // syncptracks); false in argument makes time inclusion false by default
+        bool incltime = input->getInput("include time index",sc_options,false);
         if (!scph && !scel && !scpo) {
             return 0;
         }
@@ -144,6 +160,7 @@ extern "C" {
         result->setScorePhotons(scph);
         result->setScoreElectrons(scel);
         result->setScorePositrons(scpo);
+        result->setIncludeTime(incltime); // incltime boolean is set from aquired input for the trackscoring object (sets m_include_time)
         result->setFirstEvent(first);
         result->setLastEvent(last);
         result->setBufferSize(bufSize);

HEN_HOUSE/egs++/ausgab_objects/egs_track_scoring/egs_track_scoring.h

@@ -24,6 +24,7 @@
 #  Author:          Iwan Kawrakow, 2009
 #
 #  Contributors:    Georgi Gerganov
+#                   Alexandre Demelo
 #
 ###############################################################################
 */
@@ -78,6 +79,7 @@ This ausgab object is specified via
     score photons   = yes or no # optional, yes assumed if missing
     score electrons = yes or no # optional, yes assumed if missing
     score positrons = yes or no # optional, yes assumed if missing
+    include time index = yes or no #optional, no assumed if missing
     start scoring   = event_number # optional, 0 assumed if missing
     stop  scoring   = event_number # optional, 1024 assumed if missing
     buffer size     = size         # optional, 1024 assumed if missing
@@ -116,11 +118,13 @@ class EGS_TRACK_SCORING_EXPORT EGS_TrackScoring : public EGS_AusgabObject {
             }
             else if (iarg == EGS_Application::BeforeTransport) {
                 int q = app->top_p.q;
+                EGS_Float timeindex = app->getTimeIndex(); // get current time index from the source (through the application.cpp)
                 if ((q ==  0 && m_score_photons) ||
                         (q == -1 && m_score_electrons) ||
                         (q ==  1 && m_score_positrons)) {
                     m_pts->startNewTrack(np);
                     m_pts->setCurrentParticleInfo(new EGS_ParticleTrack::ParticleInfo(q));
+                    m_pts->setCurrentTimeIndex(timeindex); // set current time index for the particle track container object
                     m_pts->addVertex(np,new EGS_ParticleTrack::Vertex(app->top_p.x,app->top_p.E));
                 }
             }
@@ -158,6 +162,9 @@ class EGS_TRACK_SCORING_EXPORT EGS_TrackScoring : public EGS_AusgabObject {
     void setScorePositrons(bool score) {
         m_score_positrons = score;
     };
+    void setIncludeTime(bool score)   {
+        m_include_time = score;
+    };
     void setFirstEvent(EGS_I64 first)  {
         m_start = first;
     };
@@ -188,6 +195,7 @@ class EGS_TRACK_SCORING_EXPORT EGS_TrackScoring : public EGS_AusgabObject {
     bool                        m_score_photons;      //!< Score photon tracks?
     bool                        m_score_electrons;    //!< Score electron tracks?
     bool                        m_score_positrons;    //!< Score positron tracks?
+    bool                        m_include_time;       //!< include time index in tracks file?
 
     string                      m_fnExtra;            //!< String to append to output file name
 

HEN_HOUSE/egs++/egs_application.cpp

@@ -28,6 +28,7 @@
 #                   Blake Walters
 #                   Reid Townson
 #                   Hubert Ho
+#                   Alexandre Demelo
 #
 ###############################################################################
 */
@@ -47,6 +48,7 @@
 #include "egs_base_source.h"
 #include "egs_simple_container.h"
 #include "egs_ausgab_object.h"
+#include "egs_base_geometry.h"
 
 #include <cstring>
 #include <cstdio>
@@ -925,8 +927,14 @@ int EGS_Application::simulateSingleShower() {
                        " attempts\n");
             return 1;
         }
+        setTimeIndex(-1);
         current_case =
             source->getNextParticle(rndm,p.q,p.latch,p.E,p.wt,p.x,p.u);
+
+        // For dynamic geometries, update positions according to the current
+        // time index, which may have been set by getNextParticle
+        geometry->getNextGeom(rndm);
+
         ireg = geometry->isWhere(p.x);
         if (ireg < 0) {
             EGS_Float t = veryFar;

HEN_HOUSE/egs++/egs_application.h

@@ -27,6 +27,7 @@
 #                   Ernesto Mainegra-Hing
 #                   Blake Walters
 #                   Reid Townson
+#                   Alexandre Demelo
 #
 ###############################################################################
 */
@@ -694,16 +695,23 @@ class EGS_EXPORT EGS_Application {
         geometry->getLabelRegions(str, regs);
     }
 
-    /*! \brief Returns the value of the \a mu synchronization parameter
+    /*! \brief Returns the value of the \a time synchronization parameter
 
-      The parameter, \a mu, is a random number on \a [0,1) associated with each
+      The parameter, \a time, is a random number on \a [0,1) associated with each
       primary history and is retrieved from \a source.  It can be used to
-      synchronize geometric parameters throughout a simulation.  If \a mu is
-      not available in \a source (i.e., the \a getMu function has not been
+      synchronize geometric parameters throughout a simulation.  If \a time is
+      not available in \a source (i.e., the \a getTime function has not been
       reimplemented in \a source), then this returns -1.
      */
-    EGS_Float getMU() {
-        return source->getMu();
+    EGS_Float getTimeIndex() {
+        return source->getTimeIndex();
+    }
+
+    /*! Sets the value of the time synchronization parameter. This will mainly
+     * be used by the dynamic geometry if it does not receive time from a source
+     * */
+    void setTimeIndex(EGS_Float temp_time) {
+        source->setTimeIndex(temp_time);
     }
 
     /*! \brief User scoring function for accumulation of results and VRT implementation

HEN_HOUSE/egs++/egs_base_geometry.h

@@ -29,6 +29,7 @@
 #                   Reid Townson
 #                   Ernesto Mainegra-Hing
 #                   Hugo Bouchard
+#                   Alexandre Demelo
 #
 ###############################################################################
 */
@@ -43,6 +44,7 @@
 #define EGS_BASE_GEOMETRY_
 
 #include "egs_vector.h"
+#include "egs_rndm.h"
 
 #include <string>
 #include <vector>
@@ -149,6 +151,21 @@ class EGS_EXPORT EGS_BaseGeometry {
     */
     virtual int isWhere(const EGS_Vector &x) = 0;
 
+    /* getNextGeom is the equivalent of getNextParticle but for the simulation
+     * object. Its goal is to determine the next state of the geometry, either
+     * by synchronizing itself to the source time parameter, or by sampling its
+     * own time parameter and updating itself accordingly if the source has
+     * provided no time index.
+     *
+     * This function has a non-empty implementation in 2 cases.
+     *
+     * 1) it is reimplemented in any composite geometry, where it will call next
+     *    geom on all of its components
+     *
+     * 2) it is reimplemented in the dynamic geometry class. This is where the
+     *    code will find the current (non static) state of the geometry. */
+    virtual void getNextGeom(EGS_RandomGenerator *rndm) {};
+
     /*! \brief Find the bin to which \a xp belongs, given \a np bin edges \a p
 
         This static method uses a binary search and is provided for the
@@ -679,7 +696,6 @@ class EGS_EXPORT EGS_BaseGeometry {
         return ++nref;
     };
 
-
     /*! \brief Decrease the reference count to this geometry
 
       Composite geometries should use this method to decrease
@@ -731,6 +747,18 @@ class EGS_EXPORT EGS_BaseGeometry {
     /*! \brief Set the labels from an input string */
     int setLabels(const string &inp);
 
+    virtual void updatePosition(EGS_Float time) { };
+
+    /* This method is essentially used to determine whether the simulation
+     * geometry contains a dynamic geometry. Like getNextGeom(), the only
+     * non-empty implementations of this function are in composite geometries
+     * (where it simply calls containsDynamic on its components), and in the
+     * dynamic geometry, where it will update the boolean reference to true and
+     * call on its base geometry. This function was conceived to be used in the
+     * view/viewcontrol (to determine whether time index objects are visible or
+     * hidden) */
+    virtual void containsDynamic(bool &hasdynamic) { };
+
 protected:
 
     /*! \brief Number of local regions in this geometry
@@ -763,7 +791,6 @@ class EGS_EXPORT EGS_BaseGeometry {
      */
     int med;
 
-
     /*! \brief Does this geometry have relative mass density scvaling?
 
      */
@@ -1116,15 +1143,13 @@ struct EGS_GeometryIntersections {
     \until make_depend
     That's all.
 
-
     Here is the complete source code of the EGS_Box class.<br>
     The header file:
     \include geometry/egs_box/egs_box.h
     The .cpp file:
     \include geometry/egs_box/egs_box.cpp
     The Makefile:
     \include geometry/egs_box/Makefile
-
 */
 
 /* \example geometry/example1/geometry_example1.cpp
@@ -1336,5 +1361,4 @@ struct EGS_GeometryIntersections {
   and the complete implementation:
 */
 
-
 #endif