Fail to histogram data after coord transform using a custom-made transformation graph

[koshchii]

I get a "DimensionError: View over subspan can only be created for contiguous range of data." when I am trying to create a 2D histogram of data, obtained upon performing coord_transform using a self-made transformation graph.

The minimal code below can be used to reproduce the error:

import scipp as sc
import scippneutron as scn
import numpy as np

def Q_elements_from_wavelength_inelastic(wavelength_i, wavelength_f, incident_beam, scattered_beam):
    e_i = incident_beam / sc.norm(incident_beam)
    e_f = scattered_beam / sc.norm(scattered_beam)
    k_i = 2 * np.pi / wavelength_i * e_i
    k_f = 2 * np.pi / wavelength_f * e_f
    q_vec = k_i - k_f
    return {'Qx': q_vec.fields.x, 'Qy': q_vec.fields.y, 'Qz': q_vec.fields.z}

Q_proj_graph = {
    ('Qx', 'Qy', 'Qz'): Q_elements_from_wavelength_inelastic
}

num_spectra = 10
lam_i = sc.scalar(1.3, unit='A')
lam_f = sc.Variable(dims=['spectrum'], unit='A', values=np.random.rand(num_spectra))
inc_beam = sc.vector(value=[6.05, 0, 0], unit='m')
scat_beam = sc.vectors(dims=['spectrum'], values=np.random.rand(num_spectra, 3), unit='m')
intensity_test = sc.Variable(dims=['spectrum'], unit='arb. units', values=np.random.rand(num_spectra))

data_test = sc.DataArray(data=intensity_test,
                         coords={'wavelength_f': lam_f,
                                 'wavelength_i': lam_i,
                                 'incident_beam': inc_beam,
                                 'scattered_beam': scat_beam,
                                })

transformed_table_test = data_test.transform_coords(["Qx","Qy","Qz"], graph=Q_proj_graph)

transformed_table_test.hist(Qx=2, Qy=2)

When I import Q_elements_from_wavelength method from scippneutron

from scippneutron.conversion.tof import Q_elements_from_wavelength

and use it instead of the self-defined Q_elements_from_wavelength_inelastic method (data_test DataArray needs to be updated correspondingly to keep only one 'wavelength' coordinate), then the histogramming works as expected.

[SimonHeybrock]

Hi @koshchii!
Thanks for reporting. The problem is that properties such as fields.x return a view, and some part of the subsequent code does not like the non-contiguous (stride 3) data. You can work around this by making copies of the fields:

    return {'Qx': q_vec.fields.x.copy(), 'Qy': q_vec.fields.y.copy(), 'Qz': q_vec.fields.z.copy()}

or an equivalent reformulation of the momentum-vector math into component-based math.

We have fixed a similar problem in the past in scipp/scipp#3255, but it seems we missed another case. I have opened an issue upstream at scipp/scipp#3660.

[koshchii]

Hi @SimonHeybrock, thanks much for getting back! The work-around solution works well ☺️

[SimonHeybrock]

Closing in favor of the upstream issue then, please subscribe there if you would like update notifications.