BUG: calculate default depth for off axis projections #5081

chrishavlin · 2024-12-12T21:23:01Z

Pretty sure this won't break the SPH rendering.... we'll see how the tests go. but would be good to have confirmation from @nastasha-w

chrishavlin · 2024-12-12T21:24:02Z

yt/visualization/plot_window.py

@@ -101,15 +105,7 @@ def get_oblique_window_parameters(

    w = tuple(el.in_units("code_length") for el in width)
    bounds = tuple(((2 * (i % 2)) - 1) * w[i // 2] / 2 for i in range(len(w) * 2))
-    if get3bounds and depth is None:


if we initially calculate a depth equal to the diagonal + a bit, we don't need this extra logic here as the off-axis projection calls to this function will always have a depth.

chrishavlin · 2024-12-12T21:25:04Z

yt/visualization/fixed_resolution.py

@@ -639,7 +639,6 @@ def _generate_image_and_mask(self, item) -> None:
                self.bounds[5] - self.bounds[4],
            )
        )
-        depth = dd.depth[0] if dd.depth is not None else None


removing this because the depth sanitization before dd creation takes care of this now.

yt/visualization/plot_window.py

matthewturk · 2024-12-13T17:28:15Z

yt/visualization/plot_window.py

+        if depth is None:
+            # off-axis projection, depth not specified
+            # -> set 'large enough' depth using half the box diagonal + margin
+            depth = np.linalg.norm(ds.domain_width.in_units("code_length")) * 1.02


Why 1.02? Does it need to be an epsilon that large?

Why 1.02?

Good question -- I was just preserving the 2% epsilon in this line

yt/yt/visualization/plot_window.py

Line 111 in e6e5f1f

bounds = bounds + (-0.51 * diag, 0.51 * diag)

Does it need to be an epsilon that large?

I suspect it could be smaller, and for periodic datasets in particular a 2% epsilon might not be ideal... @nastasha-w was there reasoning behind the magnitude of margin you used?

At least for SPH data, periodicities are applied before rotation, and are never taken to be periodic in the rotated coordinates. Therefore, the large margin shouldn't cause trouble, but a smaller value ought to be fine as well.

I guess I also figured that the half box diagonal would be a very large margin for most line of sight directions anyway

nastasha-w · 2024-12-17T13:53:37Z

This seems to make sense overall; setting a default depth somewhere should be enough, although TBH it's probably not a great idea in most situations to project a whole box off-axis; in practice that means the depth will vary across the image. It's been a while since I looked into the code and checked what's called when here, but this does seem like something that would break the tests for SPH if done wrong. (I do have a few in there that test if at least roughly the right depth is used.)

nastasha-w · 2024-12-17T14:10:53Z

Thanks for looking into this! Should we maybe add an SPH and a grid code off-axis image test to look out for this sort of thing going forward?

chrishavlin · 2024-12-17T16:26:22Z

although TBH it's probably not a great idea in most situations to project a whole box off-axis; in practice that means the depth will vary across the image.

i don't disagree necessarily, but that's what was happening implicitly here:

yt/yt/visualization/plot_window.py

Lines 107 to 111 in e6e5f1f

    
           d2 = ds.domain_width[0].in_units("code_length") ** 2 
        
           d2 += ds.domain_width[1].in_units("code_length") ** 2 
        
           d2 += ds.domain_width[2].in_units("code_length") ** 2 
        
           diag = np.sqrt(d2) 
        
           bounds = bounds + (-0.51 * diag, 0.51 * diag)

So IMO it's clearer to set that default depth explicitly on initialization of the plot window (whether the default depth is the diagonal of the whole domain or something smaller).

At least for SPH data, periodicities are applied before rotation, and are never taken to be periodic in the rotated coordinates. Therefore, the large margin shouldn't cause trouble, but a smaller value ought to be fine as well.

I'm not positive on if that's the case for gridded datasets, I'll take a closer look to be sure.

Should we maybe add an SPH and a grid code off-axis image test to look out for this sort of thing going forward?

ya, definitely -- I'm very surprised there wasn't already an off-axis image test for gridded data to catch this initially...

nastasha-w · 2024-12-17T17:11:44Z

So IMO it's clearer to set that default depth explicitly on initialization of the plot window (whether the default depth is the diagonal of the whole domain or something smaller).

Makes sense! I mostly used some of those more complicated things because I wanted to avoid breaking stuff for e.g., grid code projections. I very much appreciate a more logical/simplified approach.

ya, definitely -- I'm very surprised there wasn't already an off-axis image test for gridded data to catch this initially...

There aren't for SPH either -- I have some off-axis tests for those, but those are for fake datasets with a small number of particles. I did do a quick visual test though; I'm pasting the code below. It's for off-axis projections and slices through an SPH dataset. I'm guessing something similar would work for other hydro methods.
(I think the standard image tests were meant more as tests of the different frontends than of the different visualizations.)

# looks good!
def visual_check_offaxis_projection():
    '''
    Plot projections along the x, z axes using the axis-aligned
    methods, then check against just-off x, z and intermediate
    lines of sight to see if we get a nice, continuous rotation.
    Note that paths here are system-specific.
    '''
    simpath = '/Users/nastasha/code/ytdev_misc/sample_data/gizmo_mhd_mwdisk/'
    simfile = 'gizmo_mhd_mwdisk.hdf5'
    outdir = '/Users/nastasha/code/ytdev_misc/testimgs/'
    outfilen_temp = 'sph_projection_test_lineofsight_{projdir}.pdf'

    nsteps = 9
    projdirs = ['x', 'z'] + [f'step{i}' for i in range(nsteps)]
    epsilon = 1e-4
    north_vector = [0., 1., 0.] #y-axis
    
    ds = yt.load(simpath + simfile)
    source = ds.all_data()
    center = [0., 0., 0.]

    for projdir in projdirs:
        if projdir in ['x', 'y', 'z']:
            outfilen = outdir + outfilen_temp.format(projdir=projdir)
            projaxis = projdir
            title = f'line of sight along the {projdir}-axis'
            #dle = ds.domain_left_edge.to("kpc")
            center = ds.arr(center, 'kpc')
            le = center - ds.quan(30., 'kpc')
            re = center + ds.quan(30, 'kpc')
            reg = YTRegion(center.to('code_length'), le.to('code_length'),
                           re.to('code_length'), ds=ds)
            # need to specify width *and* region because some giant SPH
            # particle intersects the central region
            prj = yt.ProjectionPlot(ds, projaxis, ("gas", "density"),
                                    weight_field=None,
                                    buff_size=(400, 400),
                                    center=center,
                                    width=(60, "kpc"),
                                    data_source=reg)
        else:
            _projdir = projdir + f'of{nsteps}'
            outfilen = outdir + outfilen_temp.format(projdir=_projdir)
            stepnum = int(projdir[4:])
            angle = stepnum / (nsteps - 1.) * 0.5 * np.pi
            projaxis = [np.cos(angle + epsilon), 0., np.sin(angle + epsilon)]
            title = f'line of sight x- to z-axis: step {stepnum} / {nsteps}'
            prj = yt.ProjectionPlot(ds, projaxis, ("gas", "density"),
                                    width=(60, "kpc"),
                                    weight_field=None,
                                    buff_size=(400, 400),
                                    center=center,
                                    data_source=source,
                                    north_vector=north_vector,
                                    depth=(60, 'kpc'))
        
        prj.annotate_title(title)
        prj.save(outfilen)

# looks ok, seems consisitent with projections
def visual_check_offaxis_slice():
    '''
    Plot projections along the x, z axes using the axis-aligned
    methods, then check against just-off x, z and intermediate
    lines of sight to see if we get a nice, continuous rotation.
    Note that paths here are system-specific.
    '''
    simpath = '/Users/nastasha/code/ytdev_misc/sample_data/gizmo_mhd_mwdisk/'
    simfile = 'gizmo_mhd_mwdisk.hdf5'
    outdir = '/Users/nastasha/code/ytdev_misc/testimgs/'
    outfilen_temp = 'sph_slice_test_lineofsight_{projdir}.pdf'
    # bounding box recommended at loading
    simbounds = [[-335.41967773, -305.92471313, -304.52108765],
                 [283.54318237, 300.88858032, 379.63751221]]

    nsteps = 9
    projdirs = ['x', 'z'] + [f'step{i}' for i in range(nsteps)]
    epsilon = 1e-4
    north_vector = [0., 1., 0.] #y-axis
    
    ds = yt.load(simpath + simfile, bounding_box=simbounds)
    center = ([0., 0., 0.], 'kpc')
    width = (60., "kpc")

    for projdir in projdirs:
        if projdir in ['x', 'y', 'z']:
            outfilen = outdir + outfilen_temp.format(projdir=projdir)
            projaxis = projdir
            title = f'line of sight along the {projdir}-axis'
            slc = yt.SlicePlot(ds, projaxis, ("gas", "density"),
                               center=center, width=width,
                               buff_size=(400, 400),)
        else:
            _projdir = projdir + f'of{nsteps}'
            outfilen = outdir + outfilen_temp.format(projdir=_projdir)
            stepnum = int(projdir[4:])
            angle = stepnum / (nsteps - 1.) * 0.5 * np.pi
            projaxis = [np.cos(angle + epsilon), 0., np.sin(angle + epsilon)]
            title = f'line of sight x- to z-axis: step {stepnum} / {nsteps}'
            slc = yt.SlicePlot(ds, projaxis, ("gas", "density"),
                               center=center, width=width,
                               buff_size=(400, 400), north_vector=north_vector)
        slc.annotate_title(title)
        slc.save(outfilen)


if __name__ == '__main__':
    visual_check_offaxis_projection()
    visual_check_offaxis_slice()

chrishavlin · 2024-12-17T21:47:25Z

Thanks for the example tests! I think I'll open an issue for adding off axis projections tests so we can get this bug fix in a bit faster. I think I'd also like to address the magnitude of the margin in a follow up as well, @matthewturk how do you feel about that?

matthewturk · 2024-12-17T22:15:32Z

sounds good to me

chrishavlin · 2024-12-17T22:25:12Z

Ok, so I think that this one is good to go then, but needs an official approval from someone (@matthewturk @nastasha-w or anyone)

calculate default depth for off axis proj initially

f014741

chrishavlin added bug viz: 2D labels Dec 12, 2024

chrishavlin added this to the 4.4.1 milestone Dec 12, 2024

chrishavlin requested a review from nastasha-w December 12, 2024 21:23

chrishavlin commented Dec 12, 2024

View reviewed changes

yt/visualization/plot_window.py Outdated Show resolved Hide resolved

use np.linalg.norm

cb7035b

matthewturk reviewed Dec 13, 2024

View reviewed changes

chrishavlin mentioned this pull request Dec 17, 2024

Add tests for off-axis projections #5086

Open

Provide feedback

Saved searches

Use saved searches to filter your results more quickly

BUG: calculate default depth for off axis projections #5081

BUG: calculate default depth for off axis projections #5081

chrishavlin commented Dec 12, 2024 •

edited

Loading

chrishavlin Dec 12, 2024

chrishavlin Dec 12, 2024

matthewturk Dec 13, 2024

chrishavlin Dec 13, 2024

nastasha-w Dec 17, 2024

nastasha-w Dec 20, 2024

nastasha-w commented Dec 17, 2024

nastasha-w commented Dec 17, 2024

chrishavlin commented Dec 17, 2024

nastasha-w commented Dec 17, 2024

chrishavlin commented Dec 17, 2024

matthewturk commented Dec 17, 2024

chrishavlin commented Dec 17, 2024

BUG: calculate default depth for off axis projections #5081

Are you sure you want to change the base?

BUG: calculate default depth for off axis projections #5081

Conversation

chrishavlin commented Dec 12, 2024 • edited Loading

chrishavlin Dec 12, 2024

Choose a reason for hiding this comment

chrishavlin Dec 12, 2024

Choose a reason for hiding this comment

matthewturk Dec 13, 2024

Choose a reason for hiding this comment

chrishavlin Dec 13, 2024

Choose a reason for hiding this comment

nastasha-w Dec 17, 2024

Choose a reason for hiding this comment

nastasha-w Dec 20, 2024

Choose a reason for hiding this comment

nastasha-w commented Dec 17, 2024

nastasha-w commented Dec 17, 2024

chrishavlin commented Dec 17, 2024

nastasha-w commented Dec 17, 2024

chrishavlin commented Dec 17, 2024

matthewturk commented Dec 17, 2024

chrishavlin commented Dec 17, 2024

chrishavlin commented Dec 12, 2024 •

edited

Loading