import yt
import numpy as np
ds = yt.load('/Users/goldbaum/Documents/test/enzo_cosmology_plus/DD0046/DD0046')
ad = ds.all_data()
ad.set_field_parameter('omega_baryon', .02)
hI_plot = yt.PhasePlot(ad, 'baryon_overdensity', 'temperature', 'H_p0_density')
hI_profile = np.array(hI_plot.profile['H_p0_density'])
x_data = np.array(hI_plot.profile.x)
y_data = np.array(hI_plot.profile.y)
hI_plot
n_cells_plot = yt.PhasePlot(ad, 'baryon_overdensity', 'temperature', 'ones', weight_field=None)
n_cells_profiles = n_cells_plot.profile['ones']
n_cells_plot
n_cells_profiles.shape
(128, 128)
hI_profile.shape
(128, 128)
plot_object = hI_plot.plots['H_p0_density']
wh = np.where(n_cells_profiles < 20)
%matplotlib inline
from matplotlib import pyplot as plt
from matplotlib.colors import LogNorm
# need to transpose hI_profile to get the correct orientation, see the docstrings for pcolor
plt.pcolormesh(x_data, y_data, hI_profile.T, norm=LogNorm(), vmin=5e-42, vmax=5e-34, cmap='viridis')
plt.xscale('log')
plt.yscale('log')
hI_profile[wh] = 0
plt.pcolormesh(x_data, y_data, hI_profile.T, norm=LogNorm(), vmin=5e-42, vmax=5e-34, cmap='viridis')
plt.xscale('log')
plt.yscale('log')