I’m working on some updates to my yt-based pyXSIM package and I’m stuck and I thought I would email this list for some guidance.

I’m thinking of a situation where we have a data object (sphere, box, whatever) and it straddles a periodic boundary. I want to convert the coordinates such that the coordinates are translated with respect to some arbitrary origin (say, the center of a sphere but in theory it could be anywhere), but are also continuous, i.e, they do not wrap at the boundary.

I’ve looked at the functions get_radius and get_periodic_rvec in yt/fields/field_functions.py, and based on that I have come up with the code below, but it doesn’t quite work for any arbitrary value of the “ctr” argument (the origin). I think it’s because the functions I mentioned only need to calculate absolute value differences in order to compute a radius.

The x, y, and z arguments to the function below are the input coordinate arrays. Could anyone who is more familiar with this sort of thing point out what I should be doing?

Best,

John Z

def get_periodic_coords(ds, ctr, x, y, z):

coords = ds.arr(np.zeros((3, x.size)), "kpc")

coords[0, :] = (x - ctr[0]).to("kpc")

coords[1, :] = (y - ctr[1]).to("kpc")

coords[2, :] = (z - ctr[2]).to("kpc")

if sum(ds.periodicity) == 0:

return coords

dw = ds.domain_width.to("kpc")

for i in range(3):

if ds.periodicity[i]:

c = coords[i, :]

cdw = c - dw[i]

mins = np.argmin([np.abs(c), np.abs(cdw)], axis=0)

ii = np.where(mins == 1)[0]

coords[i, ii] = coords[i, ii] - dw[i]

return coords