"""Demonstration script to trace the evolution of one galaxy. It generates a plot of the stellar mass and metallicity of the galaxy throughout cosmic time, including its progenitors. This is an example of how to use information in the high-level tables and SpiderWeb catalogue to trace the evolution of a galaxy. It also illustrates using the read_index keyword in reading data from HDF5 files. """ # Import required packages import numpy as np import hydrangea as hy import matplotlib.pyplot as plt from astropy.cosmology import Planck13 # Set script parameters sim_index = 0 # Which simulation do we want? galaxy_id = 7679 # Which galaxy should we follow? plotloc = 'galaxy_evolution.png' # Where to save the output plot? # Prepare the plot fig = plt.figure(figsize=(5.2, 4)) # Set up the simulation object for the run we're working with, to easily # get the relevant file paths sim = hy.Simulation(index=sim_index) # ------------------------------------------------------ # I.) Simple start: the galaxy's own stellar mass growth # ------------------------------------------------------ # Load stellar masses from the high-level catalogues (in log M/M_sun). # By specifying a read_index, we only load the masses for the one specific # galaxy we are interested in (at all 30 snapshots) m_star_gal = hy.hdf5.read_data(sim.fgt_loc, 'Mstar', read_index=galaxy_id) # The galaxy may not exist in all snapshots (either because it has not yet # emerged above the resolution threshold, has been disrupted, or has been # temporarily not found by Subfind). In these cases, -1 is used as a # filler value. Let's find only snapshots with real masses: ind_good_snap = np.nonzero(m_star_gal >= 0)[0] # Plot the stellar mass growth history (against age of the Universe) snap_times = hy.snep_times(time_type='age') plt.plot(snap_times[ind_good_snap], m_star_gal[ind_good_snap], color='red', linewidth=2, marker='o', label='Main progenitor') # --------------------------------------------------------------------- # II.) Evolution of arbitrary quantities (not in high-level catalogues) # --------------------------------------------------------------------- # Now let's add a quantity that is not in the high-level catalogues # already, such as stellar metallicity. For this, we need to go back # to the individual subfind catalogues, so we need the subfind catalogue # index of the galaxy in each snapshot: shi_gal = hy.hdf5.read_data(sim.fgt_loc, 'SHI', read_index=galaxy_id) # Now load the metallicity from each individual catalogue: z_star = np.zeros(30) for isnap in range(30): # We can skip snapshots in which the galaxy does not exist if shi_gal[isnap] < 0: continue # First, get the path to the (first file of the) subfind catalogue subfind_file = sim.get_subfind_file(isnap) # Set up a reader object for the multi-file catalogue (specifying # that we only care about one particular entry, for speed-up) subhalo_galaxy = hy.SplitFile(subfind_file, 'Subhalo', read_index=shi_gal[isnap], verbose=0) # (Stellar) metallicity is stored in a separate group 'Stars', which # uses a double-underscore separator to access it as an attribute z_star[isnap] = subhalo_galaxy.Stars__Metallicity # Add the metallicity to the plot (on a second axis for clarity) ax = plt.gca() ax2 = ax.twinx() plt.axes(ax2) plt.plot(snap_times[ind_good_snap], z_star[ind_good_snap], color='royalblue', linewidth=2, linestyle='--') # --------------------------------------------------------------- # III.) Final part: total stellar mass of all progenitor galaxies # --------------------------------------------------------------- # So far, we have only considered the evolution of the galaxy itself # (i.e. its main progenitor). Let's add the evolution of the total mass # in all progenitor galaxies, including those that merge into the # target galaxy by z = 0. # Find all galaxies that merge with the target by z = 0: all those that # have the target's ID as their MergeList entry merge_ids = hy.hdf5.read_data(sim.spider_loc, 'MergeList', read_index=29, index_dim=1) ind_progenitors = np.nonzero(merge_ids == galaxy_id)[0] print(f'Galaxy {galaxy_id} has {len(ind_progenitors)} progenitors in total.') # We now need the stellar mass table again, but this time for all progenitor # galaxies (in practice, we could have loaded this directly in part I) log_m_star_progenitors = hy.hdf5.read_data(sim.fgt_loc, 'Mstar', read_index=ind_progenitors) # Sum the progenitor masses in each snapshot (taking care of filler values) m_star_progenitors = 10.0**log_m_star_progenitors ind_dummy = np.nonzero(log_m_star_progenitors < 0) m_star_progenitors[ind_dummy] = 0 m_star_tot = np.sum(m_star_progenitors, axis=0) # Plot progenitor masses on top of the galaxy's own mass evolution plt.axes(ax) plt.plot(snap_times, np.log10(m_star_tot), color='indianred', linewidth=1, label='All progenitors') # Some final embellishments on the figure plt.legend(loc=4) age_z0 = snap_times[29] ax.set_xlim((0, age_z0)) ax.set_ylim((6, 11)) ax2.set_xlim((0, age_z0)) ax2.set_ylim((0, 0.03)) ax.set_xlabel('Age of the Universe [Gyr]') ax.set_ylabel(r'$\log_{10}\,(M_\mathrm{star}\,/\,\mathrm{M}_\odot)$') ax2.set_ylabel(r'Stellar metallicity ($Z$)', labelpad=20) ax2.yaxis.label.set_rotation(-90.0) # Colour two y axes according to which data they represent ax2.spines['left'].set_color('red') ax.tick_params(axis='y', colors='red') ax.yaxis.label.set_color('red') ax2.spines['right'].set_color('royalblue') ax2.tick_params(axis='y', colors='royalblue') ax2.yaxis.label.set_color('royalblue') # Add redshifts along the upper x-axis ax3 = ax.twiny() zred = np.array((4,2,1,0.5,0.3,0.1,0)) age = Planck13.age(zred).value ax3.set_xlim((0, age_z0)) ax3.set_xticks(age) ax3.set_xticklabels([str(_z) for _z in zred]) ax3.set_xlabel('Redshift') # Save the figure plt.subplots_adjust(left=0.14, bottom=0.15, right=0.82, top=0.88) plt.show plt.savefig(plotloc, transparent=False, dpi=150)