I’m running a set of simulations to study the evolution of Super Massive black holes.

In principle, I have the black hole masses and other property called merger tree weight in Mpc^{-3} unit which is different for each merger tree created. If I want to construct the black hole mass function I should read the black holes masses, and the mergerTreeWeight but I don’t know how to implement the merger tree weight exactly in the code.

Also, I have the indexes where every merger tree start and end. I know it is useful because if the mass function consists of a set of bins in mass. Let’s say that each bin ranges from

$\log(M)-\mathrm{d}(\log M)/2$ to $\log(M)+\mathrm{d}(\log M)/2$. To compute the mass function in a bin, I’d select all black holes with mass in that range. Then, just sum the mergerTreeWeight of those galaxies, and divide by the bin width, $\mathrm{d}(\log M)$.

My problem is related to how to implement this process in a python code.

*EDIT*

I think it should be understood considering three variables such as BH = blackHoleMass , MTW= MergerTreeWeight and z= Mergertreeindex. The “z” variable give the index where every merger tree (with different weight) start `BH = [1e9, 4e7, 7e5,2e3 ] MTW= [1.2e-5, 4.6e-4,1.4e-6] z=[0,2,3] `

. BH[:3] is assigned to the first merger tree, which has a MergerTreeWeight equals to MTW[0], BH[3] is assigned to the second MergerTreeWeight, MTW[1] and BH[4] is assigned to the third MergerTreeWeight.

##
1 Answer

Assuming that the three arrays `BH`

, `MTW`

and `z`

are numpy arrays, you can do the following:

```
# First determine the length of each merger tree
MTlength = np.diff(z) # length of each MT
# The length of the last MT is missing, let's add it
MTlength = np.append(MTlength, len(BH)-z[-1])
# Now let's build an array with len(BH) elements that for each black hole
# contains the corresponding MT weight.
BHW = np.repeat(MTW, MTlength)
# Define the array of mass bins, this is an example that goes from
# 10 to 10^10 in 30 logspaced bins
mbins = np.logspace(1,10,30)
# Build the mass function
MF, _ = np.histogram(BH, bins=mbins, weights=BHW)
# Divide by the size of each bin
MF = MF/np.diff(mbins)
```