Skip to content
This repository was archived by the owner on Aug 28, 2024. It is now read-only.

Commit

Permalink
The encounter output files are now merged properly into the .enc inst…
Browse files Browse the repository at this point in the history
…ance variable
  • Loading branch information
daminton committed Dec 6, 2022
1 parent d34d36e commit dd30791
Show file tree
Hide file tree
Showing 4 changed files with 12 additions and 10 deletions.
12 changes: 7 additions & 5 deletions examples/Fragmentation/Fragmentation_Movie.py
Original file line number Diff line number Diff line change
Expand Up @@ -84,7 +84,8 @@ class AnimatedScatter(object):
"""An animated scatter plot using matplotlib.animations.FuncAnimation."""

def __init__(self, sim, animfile, title, nskip=1):
self.ds = sim.enc.mean(dim="encounter") # Reduce out the encounter dimension
tgood = sim.enc.where(sim.enc['Gmass'] > 9e-8).time
self.ds = sim.enc.sel(time=tgood)
nframes = int(self.ds['time'].size)
self.sim = sim
self.title = title
Expand All @@ -109,8 +110,8 @@ def setup_plot(self):

# Calculate the distance along the y-axis between the colliding bodies at the start of the simulation.
# This will be used to scale the axis limits on the movie.
rhy1 = self.ds['rh'].isel(time=0).sel(name="Body1",space='y').values[()]
rhy2 = self.ds['rh'].isel(time=0).sel(name="Body2",space='y').values[()]
rhy1 = self.ds['rh'].sel(name="Body1",space='y').isel(time=0).values[()]
rhy2 = self.ds['rh'].sel(name="Body2",space='y').isel(time=0).values[()]

scale_frame = abs(rhy1) + abs(rhy2)
ax = plt.Axes(fig, [0.1, 0.1, 0.8, 0.8])
Expand All @@ -133,12 +134,13 @@ def center(Gmass, x, y):
x = x[~np.isnan(x)]
y = y[~np.isnan(y)]
Gmass = Gmass[~np.isnan(Gmass)]
x = x[Gmass]
x_com = np.sum(Gmass * x) / np.sum(Gmass)
y_com = np.sum(Gmass * y) / np.sum(Gmass)
y_com = #np.sum(Gmass * y) / np.sum(Gmass)
return x_com, y_com

Gmass, rh, point_rad = next(self.data_stream(frame))
x_com, y_com = center(Gmass, rh[:,0], rh[:,1])
#x_com, y_com = center(Gmass, rh[:,0], rh[:,1])
self.scatter_artist.set_offsets(np.c_[rh[:,0] - x_com, rh[:,1] - y_com])
self.scatter_artist.set_sizes(point_rad**2)
return self.scatter_artist,
Expand Down
2 changes: 1 addition & 1 deletion python/swiftest/swiftest/simulation_class.py
Original file line number Diff line number Diff line change
Expand Up @@ -2770,7 +2770,7 @@ def _preprocess(ds, param):
return io.process_netcdf_input(ds,param)
partial_func = partial(_preprocess, param=self.param)

self.enc = xr.open_mfdataset(enc_files,parallel=True,combine="nested", concat_dim="encounter",preprocess=partial_func)
self.enc = xr.open_mfdataset(enc_files,parallel=True,preprocess=partial_func,mask_and_scale=True)
self.enc = io.process_netcdf_input(self.enc, self.param)

return
Expand Down
4 changes: 2 additions & 2 deletions src/symba/symba_io.f90
Original file line number Diff line number Diff line change
Expand Up @@ -137,7 +137,7 @@ end subroutine symba_io_encounter_initialize_output
module subroutine symba_io_encounter_write_frame(self, nc, param)
!! author: David A. Minton
!!
!! Write a frame of output of an encounter list structure.
!! Write a frame of output of an encounter trajectory.
use netcdf
implicit none
! Arguments
Expand Down Expand Up @@ -369,7 +369,7 @@ module subroutine symba_io_start_encounter(self, param, t)
call self%encounter_history%reset()

! Empty out the time slot array for the next pass
self%encounter_history%tvals(:) = -huge(1.0_DP)
self%encounter_history%tvals(:) = huge(1.0_DP)

! Take the snapshot at the start of the encounter
call self%snapshot(param, t)
Expand Down
4 changes: 2 additions & 2 deletions src/symba/symba_util.f90
Original file line number Diff line number Diff line change
Expand Up @@ -1313,7 +1313,7 @@ module subroutine symba_util_take_encounter_snapshot(self, param, t)
!! author: David A. Minton
!!
!! Takes a minimal snapshot of the state of the system during an encounter so that the trajectories
!! Can be played back through the encounter
!! can be played back through the encounter
implicit none
! Internals
class(symba_nbody_system), intent(inout) :: self !! SyMBA nbody system object
Expand Down Expand Up @@ -1401,7 +1401,7 @@ module subroutine symba_util_take_encounter_snapshot(self, param, t)
! Find out which time slot this belongs in by searching for an existing slot
! with the same value of time or the first available one
do i = 1, self%encounter_history%nframes
if (t >= self%encounter_history%tvals(i)) then
if (t <= self%encounter_history%tvals(i)) then
snapshot%tslot = i
self%encounter_history%tvals(i) = t
exit
Expand Down

0 comments on commit dd30791

Please sign in to comment.