Made minor edits to implimentation.

Moved particle track subroutines to the output module.  Particle tracks can now be written in any run mode that calls transport.  track.py should correctly handle any number of input particle tracks.

src/DEPENDENCIES

@@ -297,7 +297,6 @@ physics.o: material_header.o
 physics.o: mesh.o
 physics.o: output.o
 physics.o: particle_header.o
-physics.o: particle_restart.o
 physics.o: particle_restart_write.o
 physics.o: random_lcg.o
 physics.o: search.o

src/output.F90

@@ -1945,4 +1945,44 @@ contains
 
   end function get_label
 
+!===============================================================================
+! INITIALIZE_PARTICLE_TRACK opens a particle track output file.
+! 
+! TODO: This subroutine needs to be modified to work with HDF5 files.  It
+! should also probably write a header that identifies the file as a particle
+! track and maybe adds particle identifying information (batch #, etc.).
+!===============================================================================
+
+  subroutine initialize_particle_track()
+    character(MAX_FILE_LEN) :: filename
+
+    filename = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
+         // '_' // trim(to_str(current_work)) // '.binary'
+    open(UNIT=UNIT_TRACK, FILE=filename, STATUS='replace', &
+         ACCESS='stream')
+
+  end subroutine initialize_particle_track
+
+!===============================================================================
+! WRITE_PARTICLE_TRACK outputs particle position to a binary file.
+! 
+! TODO: This subroutine needs to be modified to work with HDF5 files.  Perhaps
+! it should also be made somehow more general so that it can output
+! information other than just particle position.
+!===============================================================================
+
+  subroutine write_particle_track()
+    write(UNIT_TRACK) p % coord0 % xyz
+  end subroutine write_particle_track
+
+!===============================================================================
+! FINALIZE_PARTICLE_TRACK closes the particle track file.
+! 
+! TODO: This subroutine needs to be modified to work with HDF5 files.
+!===============================================================================
+
+  subroutine finalize_particle_track()
+    close(UNIT=UNIT_TRACK)
+  end subroutine finalize_particle_track
+
 end module output

src/particle_restart.F90

@@ -11,7 +11,6 @@ module particle_restart
   use physics,         only: transport
   use random_lcg,      only: set_particle_seed
   use source,          only: initialize_particle
-  use string,          only: to_str
 
 #ifdef HDF5
   use hdf5_interface 
@@ -20,7 +19,6 @@ module particle_restart
   implicit none
   private
   public ::  run_particle_restart
-  public ::  write_particle_track
 
 #ifdef HDF5
   integer(HID_T) :: hdf5_particle_file
@@ -144,36 +142,12 @@ contains
          current_gen - 1)*n_particles + p % id
     call set_particle_seed(particle_seed)
 
-    ! Open particle track output file.
-    if (write_track) then
-      filename = trim(path_output) // 'track_' // trim(to_str(current_batch)) &
-           // '_' // trim(to_str(current_work)) // '.binary'
-      open(UNIT=UNIT_TRACK, FILE=filename, STATUS='replace', &
-           ACCESS='stream')
-    end if
-
     ! transport neutron
     call transport()
 
-    ! Close particle track output file.
-    if (write_track) then
-      close(UNIT=UNIT_TRACK)
-    endif
-
     ! write output if particle made it
     call print_particle()
 
   end subroutine run_particle_restart
 
-!===============================================================================
-! WRITE_PARTICLE_TRACK outputs particle position to a binary file.  This
-! subroutine needs to be modified to work with HDF5 files.  Perhaps it should
-! also be made somehow more general so that it can output information other
-! than just particle position.
-!===============================================================================
-
-  subroutine write_particle_track()
-    write(UNIT_TRACK) p % coord0 % xyz
-  end subroutine write_particle_track
-
 end module particle_restart

src/physics.F90

@@ -13,9 +13,11 @@ module physics
   use interpolation,          only: interpolate_tab1
   use material_header,        only: Material
   use mesh,                   only: get_mesh_indices
-  use output,                 only: write_message
+  use output,                 only: write_message, initialize_particle_track, &
+                                    write_particle_track, &
+                                    finalize_particle_track
   use particle_header,        only: LocalCoord
-  use particle_restart,       only: write_particle_track
+!  use particle_restart,       only: write_particle_track
   use particle_restart_write, only: write_particle_restart
   use random_lcg,             only: prn
   use search,                 only: binary_search
@@ -74,6 +76,11 @@ contains
     ! Force calculation of cross-sections by setting last energy to zero 
     micro_xs % last_E = ZERO
 
+    ! Prepare to write out particle track.
+    if (write_track) then
+      call initialize_particle_track()
+    endif
+
     do while (p % alive)
 
       ! Write particle track.
@@ -106,9 +113,11 @@ contains
         coord % xyz = coord % xyz + distance * coord % uvw
         coord => coord % next
       end do
+
+      ! Write particle track
       
-      if (run_mode == MODE_PARTICLE) then
-        write(UNIT_TRACK) p % coord0 % xyz
+      if (write_track) then
+        call write_particle_track()
       end if
 
       ! Score track-length tallies
@@ -185,6 +194,11 @@ contains
 
     end do
 
+    ! Finish particle track.
+    if (write_track) then
+      call finalize_particle_track()
+    endif
+
   end subroutine transport
 
 !===============================================================================

src/utils/track.py

@@ -1,7 +1,19 @@
 #!/usr/bin/env python2
 """Convert binary particle track to VTK poly data.
 
-Run 'track.py -help' for usage.
+Usage information can be obtained by running 'track.py --help':
+
+    usage: track.py [-h] [-o OUT] IN [IN ...]
+
+    Convert binary particle track file to a .pvtp file.
+
+    positional arguments:
+      IN                    Input particle track data filename(s).
+
+    optional arguments:
+      -h, --help            show this help message and exit
+      -o OUT, -out OUT, --out OUT
+                            Output VTK poly data filename.
 
 """
 
@@ -17,13 +29,14 @@ def _parse_args():
         description='Convert binary particle track file to a .pvtp file.')
     parser.add_argument('input', metavar='IN', type=argparse.FileType('rb'),
                         nargs='+',
-                        help='Input particle track data filename.')
-    parser.add_argument('-o', '-out', '--out', metavar='OUT', type=str, dest='out',
+                        help='Input particle track data filename(s).')
+    parser.add_argument('-o', '-out', '--out', metavar='OUT', type=str,
+                        dest='out',
                         help='Output VTK poly data filename.')
 
     # Parse and return commandline arguments.
     return parser.parse_args()
-    
+
 
 def main():
     # Parse commandline argumetns.
@@ -31,7 +44,7 @@ def main():
     
     # Make sure that the output filename ends with '.pvtp'.
     if not args.out:
-        args.out = ''.join([os.path.splitext(args.input.name)[0], '.pvtp'])
+        args.out = 'tracks.pvtp'
     elif os.path.splitext(args.out)[1] != '.pvtp':
         args.out = ''.join([args.out, '.pvtp'])
 
@@ -39,18 +52,14 @@ def main():
     cells = vtk.vtkCellArray()
     j = 0
     k = 0
-    arr = vtk.vtkIntArray()
     for fin in args.input:
         track = fin.read()
         coords = [struct.unpack("ddd", track[24*i : 24*(i+1)])
                   for i in range(len(track)/24)]
         n_points = len(coords)
-
-        i=0
+        
         for triplet in coords:
             points.InsertNextPoint(triplet)
-            arr.InsertTuple1(i,k)
-            i+=1
         
         # Create VTK line and assign points to line.
         line = vtk.vtkPolyLine()
@@ -61,8 +70,7 @@ def main():
         cells.InsertNextCell(line)
         j += n_points
         k += 1
-
-    points.SetData(arr)
+    global data
     data = vtk.vtkPolyData()
     data.SetPoints(points)
     data.SetLines(cells)
@@ -78,9 +86,7 @@ if __name__ == '__main__':
     main()
 
 
-#### The following code is retained for debugging purposes.  If 'data' in the
-#### in the function 'main' is made a global variable, the following code will
-#### automatically plot the tracks.
+#### The following code is retained for debugging purposes:
 
 ##poly_mapper = vtk.vtkPolyDataMapper()
 ##poly_mapper.SetInputConnection(data.GetProducerPort())