yt-project/yt

unexpected output when using ParticlePhasePlot

Horizonatpku opened this issue · 2 comments

Bug report

Bug summary
hi everyone,

I am trying to use the function ParticlePhasePlot to plot particles' quantities. However, some unexpected results come out.

I tried to plot (x-R-"d2center" , y-\Phi-"d2phi" , z-R-"d2center"), and the particles are colored by z-R-"d2center".
Ideally, the color should vary monotonically along the x-axis, since the z - R is equal to x - R
But the output figure shows that there are particles colored with red dispersed in the full figure.

Actual outcome

testPTC_2d-Profile_d2center_d2phi_d2center

Here is the code:
Code for reproduction

import yt
from yt.frontends import boxlib
from yt.frontends.boxlib.data_structures import AMReXDataset
import numpy as np

rootdir="./Output"
ds = AMReXDataset(rootdir+'/Plt00000')

ptc_xname="d2center"
ptc_yname="d2phi"
ptc_zname="d2center"
ptc_cmapname = "hsv"
ptcslc = yt.ParticlePhasePlot(
    ds,
    ptc_xname,
    ptc_yname,
    [ptc_zname],
    x_bins=800,
    y_bins=800,
)
ptcslc.set_log(ptc_xname,False)
ptcslc.set_log(ptc_yname,False)
ptcslc.set_log(ptc_zname,False)
ptcslc.set_cmap(("particles", ptc_zname), ptc_cmapname)
ptcslc.set_zlim(("particles", ptc_zname), 0.04, 0.065)
ptcslc.save("testPTC")

I do not know why this happens and could you please give me some suggestions, thanks a lot!

best wishes,
Yifeng He

Version Information

  • Operating System: Windows 10
  • Python Version: 3.9.13
  • yt version: 4.3.0
  • Other Libraries (if applicable): numpy

Might be related to the weight parameter -- since it's not set in your code, it defaults to None, and so values in a given bin in the phase plot will be summed. You could quickly check this by decreasing x_bins and y_bins: as you decrease the number of bins you should have more particles falling within a given bin and so more pixels maxing out at the upper end of your colorscale.

If you find that is problem, you could try weighting the phase plot by the number of particles in each bin using a "ones" field (where the field value = 1 for every particle). You might already have a field like this in your ds.derived_field_list, but if not you can add "ones" field to your dataset with

import numpy as np
def _particle_ones(field, data):
    return np.ones(data['particles', 'd2center'].shape)

ds.add_field(
    name=('particles', 'ones'),
    function=_particle_ones,
    sampling_type="particle",
    units="",
)

then use weight=('particles', 'ones') in your ParticlePhasePlot.

Let me know if that helps.

It works!
Thanks a lot!