Hi Andrew,

I think I see what you're seeing in the plots but it's very difficult to say what's happening without an example that one of us can run locally. Can you trigger this behavior using one of the public datasets from yt-project.org/data? If so I'd encourage you to file an issue on yt's github so we don't lose track of this. You can also share the dataset you're using if that's kosher with your collaborators.

You may also want to examine the image arrays and compare the numerical values in the two plots without having to squint at the colorbar, you can get an image from a ProjectionPlot like so:

image = prj.frb['gas', 'density']

image will then be a 2D array that you can query for the pixel-by-pixel values in the ProjectionPlot.

On Wed, Jul 10, 2019 at 4:06 PM Andrew Wu <andrewmw@princeton.edu> wrote:

Hi yt,

When making projection plots of total gas density and cold gas density, averaged over the projection length (i.e. with a uniform weighting), and using cut_region() to filter out hot gas to make the cold gas plots, it appears that some areas in the cold gas plot have a higher density than the corresponding areas in the total gas plot, which shouldn't happen. If you look at the top-middle of the two plots (link below) and compare the densities, the cold gas appears to have a higher density there. The code is below:

(link to plots: https://drive.google.com/drive/folders/1UJ6fmKNFQzGr8cfjov5-GwOzfS8nt2e8?usp=sharing)


import yt

ds = yt.load("/tigress/cen/CZ15/redshift0100")

z = 0.62

# the box size is 120/0.7/(1+z) Mpc, where z is redshift
boxSize = 120e3/0.7/(1+z)

# virial radius in code units, r_v is virial radius in kpc
codeWidth = r_v / boxSize

center = [xCoord, yCoord, zCoord]

left_edge = [xCoord-2.1*codeWidth, yCoord-2.1*codeWidth, zCoord-2.1*codeWidth]
right_edge = [xCoord+2.1*codeWidth, yCoord+2.1*codeWidth, zCoord+2.1*codeWidth]

# want a box with full width about 4 times the virial radius
reg = ds.region(center, left_edge,right_edge)


# total gas plot

prj = yt.ProjectionPlot(ds, "x", "density", center=center, width=(4*r_v, 'kpc'),
                        data_source=reg, method='integrate', weight_field='ones')

u = 100 / 0.048
l = 0.07 / 0.048

prj.set_zlim('density', l,u)
prj.set_unit('density','mbd')

prj.set_cmap(field='density', cmap='CMRmap')
prj.show()

# cold gas plot

cold_reg = reg.cut_region(["obj['temperature'] < 1e5"])

prj = yt.ProjectionPlot(ds, "x", "density", center=center, width=(4*r_v, 'kpc'),
                        data_source=cold_reg, method='integrate', weight_field='ones')

prj.set_zlim('density', l,u)
prj.set_unit('density','mbd')

prj.set_cmap(field='density', cmap='CMRmap')

prj.show()


Thanks,
Andrew
_______________________________________________
yt-users mailing list -- yt-users@python.org
To unsubscribe send an email to yt-users-leave@python.org