How does the Visualize Galaxies and Halos get the velocity and velocity dispersion map?

Dingyi Zhao
  • 15 Oct '23

Hi,

I'm trying to generate some velocity and velocity dispersion map of galaxies in TNG-50. However, I find that the map generated by myself are different from the maps generated by Visualize Galaxies and Halos. I wonder if there are some documents about Visualize Galaxies and Halos.

Dingyi Zhao
  • 3
  • 15 Oct '23

Here is an example of the velocity dispersion map generated by me with subhalo_id=12 and z=0.05.

image001.png

It is generated by the following python code:

pi = np.loadtxt('/Users/zhaodingyi/Downloads/p.txt').T
#pi[0] is x, pi[1] is y, pi[2] is v, pi[3] is mass for particle in subhalo 12
bins=100
rmax=20
pixel = 2*rmax/bins
def get_aver(x, z, bins, v):
    '''
    get the average value of v in the bin where (x,z) is in
    '''
    index_x = np.array((x + rmax) / pixel, dtype=int)
    index_z = np.array((z + rmax) / pixel, dtype=int)
    index_x[index_x>bins - 1] = bins - 1
    index_x[index_x<0] = 0
    index_z[index_z>bins - 1] = bins - 1
    index_z[index_z<0] = 0
    return v[index_x,index_z]
mass, X, Y = np.histogram2d(pi[0], pi[1], bins=bins, range=[[-rmax,rmax],[-rmax,rmax]], weights = pi[3])
v, X, Y = np.histogram2d(pi[0], pi[1], bins=bins, range=[[-rmax,rmax],[-rmax,rmax]], weights = pi[2]*pi[3])
v=v/mass
sigma, X, Y = np.histogram2d(pi[0], pi[1], bins=bins, range=[[-rmax,rmax],[-rmax,rmax]], weights = (pi[2] - get_aver(pi[0], pi[1], bins, v))**2*pi[3])
sigma = np.sqrt((sigma / mass))

The map shown above is very different from the velocity dispersion map generated by the TNG website, Visualize Galaxies and Halos, shown as below.

image002.png

Meanwhile, I can generate a map similar to Visualize Galaxies and Halos, with the code:

image003.png

sigma, X, Y = np.histogram2d(pi[0], pi[1], bins=bins, range=[[-rmax,rmax],[-rmax,rmax]], weights = (pi[2] - get_aver(pi[0], pi[1], bins, v.T))**2*pi[3]) ## the only difference is v.T, while in origin code is v

I also generated the velocity map by another method, it is also different from the map generated by Visualize Galaxies and Halos, but the same as mine. The code is shown below:

image004.png

def window_bin_data_2D(x,y,z,weights,f,x_step=0.02,y_step=0.02,x_window_length=0.2,y_window_length=0.2,threshold=1,xrange=None,yrange=None):
    '''
    2D sliding window function, when calculating the velocity dispersion, the x and y is the coordinate of the star particle, z is the velocity, weights is the mass, f is the function to calculate the weighted standard deviation.
    '''
    if(xrange==None):
        x_min = np.min(x)
        x_max = np.max(x)
    else:
        x_min = xrange[0]
        x_max = xrange[1]
    if(yrange==None):
        y_min = np.min(y)
        y_max = np.max(y)
    else:
        y_min = yrange[0]
        y_max = yrange[1]
    x_bin_boundary_min = np.arange(x_min,x_max,x_step)
    x_bin_boundary_max = np.arange(x_min+x_window_length,x_max+x_window_length,x_step)
    y_bin_boundary_min = np.arange(y_min,y_max,y_step)
    y_bin_boundary_max = np.arange(y_min+y_window_length,y_max+y_window_length,y_step) ## the boundaries of all the windows we need to calculate
    N_X = np.min((len(x_bin_boundary_max),len(x_bin_boundary_min)))
    N_Y = np.min((len(y_bin_boundary_max),len(y_bin_boundary_min)))
    X = []
    Y = []
    Z = []
    for i in range(N_X):
        temp_index = (x>=x_bin_boundary_min[i]) & (x<x_bin_boundary_max[i])
        temp_x = (x_bin_boundary_min[i]+x_bin_boundary_max[i])/2
        for j in range(N_Y):
            index = temp_index & (y>=y_bin_boundary_min[j]) & (y<y_bin_boundary_max[j]) ## index of the particles in the window
            X.append(temp_x)
            Y.append((y_bin_boundary_min[j]+y_bin_boundary_max[j])/2)
            if(np.sum(index)<threshold):
                Z.append(np.nan)
                continue
            Z.append(f(z[index],weights[index])) ## calculate the velocity dispersion of the particles in the window
    X = np.array(X)
    Y = np.array(Y)
    Z = np.array(Z)
    X = X.reshape(N_X,N_Y)
    Y = Y.reshape(N_X,N_Y)
    Z = Z.reshape(N_X,N_Y). ## Z is the velocity dispersion at X,Y
    return X,Y,Z

So, I wonder if Visualize Galaxies and Halos make a mistake, or if there are some bugs in my code. Would you like to give me some advice or some documents about Visualize Galaxies and Halos?

Dylan Nelson
  • 24 Oct '23

Dear Dingyi,

Thank you for the report, you are absolutely correct. In fact I realized this issue (affecting only the velsigma_los field) awhile ago, but did not update the code running on the website. I've updated it now, so you should see the results you expect for any future visualizations.

Unfortunately the analysis code which powers the "Visualize Galaxies and Halos" tool is not yet publicly available, though hopefully it will be soon. In the meantime, let me know if you have any other questions about it.

  • Page 1 of 1