Mantid SANS Qx Qy Data test

q2d.py

import matplotlib.pyplot as plt
import numpy as np

# binnned values
qx = np.linspace(0.001, 0.9, 193)
qy = np.linspace(0.001, 0.7, 257)

dqx = np.sqrt(qx)
dqy = np.sqrt(qy)

# Y cannot be binned in Mantid!!!
qy_bin_centers = (qy[1:]+qy[:-1])/2.
dqy_bin_centers = (dqy[1:]+dqy[:-1])/2.

#
qx_grid, qy_bin_centers_grid = np.meshgrid(qx, qy_bin_centers)
dqx_grid, dqy_bin_centers_grid = np.meshgrid(dqx, dqy_bin_centers)

i = np.sqrt((qx_grid-0.45)**2+(qy_bin_centers_grid-0.35)**2) * 100
i_sigma = np.sqrt(i)

# Q WS
ws_q = CreateWorkspace(
    DataX=qx_grid,  # 2D
    DataY=i,  # 2D
    DataE=i_sigma,  # 2D
    NSpec=256,
    UnitX='MomentumTransfer',
    VerticalAxisUnit='MomentumTransfer',
    VerticalAxisValues=qy_bin_centers,  # 1D
)

# dQ WS
ws_dq = CreateWorkspace(
    DataX=dqx_grid,  # 2D
    DataY=i,  # 2D
    DataE=i_sigma,  # 2D
    NSpec=256,
    UnitX='MomentumTransfer',
    VerticalAxisUnit='MomentumTransfer',
    VerticalAxisValues=dqy_bin_centers,  # 1D
)

GroupWorkspaces(InputWorkspaces=[ws_q, ws_dq], OutputWorkspace="qxqy")

# Andrei solution
# construct three MDWorkspaces: I(Qx, Qy),  dQx(Qx, Qy), and dQy(Qx, Qy).
# Then everytime we rebin I(Qx, Qy) we can rebin dQx and dQy in the same way in order to recalculate the resolutions


# binnned values
qx = np.linspace(0.001, 0.9, 192)
qy = np.linspace(0.001, 0.7, 256)

dqx = np.sqrt(qx)
dqy = np.sqrt(qy)

#
qx_grid, qy_grid = np.meshgrid(qx, qy)
dqx_grid, dqy_grid = np.meshgrid(dqx, dqy)

i = np.sqrt((qx_grid-0.45)**2+(qy_grid-0.35)**2) * 100
i_sigma = np.sqrt(i)

iqxqy = CreateMDHistoWorkspace(
    SignalInput=i.ravel(),
    ErrorInput=i_sigma.ravel(),
    Dimensionality='2',
    Extents='-1,1,-1,1',
    NumberOfBins='192,256',
    Names='Qx,Qy',
    Units='A-1,A-1'
)

dqxqxqy = CreateMDHistoWorkspace(
    SignalInput=qx_grid.ravel(),
    ErrorInput=dqx_grid.ravel(),
    Dimensionality='2',
    Extents='-1,1,-1,1',
    NumberOfBins='192,256',
    Names='Qx,Qy',
    Units='A-1,A-1'
)

dqyqxqy = CreateMDHistoWorkspace(
    SignalInput=qy_grid.ravel(),
    ErrorInput=dqy_grid.ravel(),
    Dimensionality='2',
    Extents='-1,1,-1,1',
    NumberOfBins='192,256',
    Names='Qx,Qy',
    Units='A-1,A-1'
)

GroupWorkspaces(InputWorkspaces=[iqxqy, dqxqxqy,
                                 dqyqxqy], OutputWorkspace="qxqymd")

###############################################################################
# MDHistoWorkspace canot be binned :( Do Andrei solution but with Workspace2D
#  I(Qx, Qy),  dQx(Qx, Qy), and dQy(Qx, Qy)

# binned values
qx_bin = np.linspace(0.001, 0.9, 193)
qy_bin = np.linspace(0.001, 0.7, 257)

dqx_bin = np.sqrt(qx_bin)
dqy_bin = np.sqrt(qy_bin)

# Binned values
qx_bin_centers = (qx_bin[1:]+qx_bin[:-1])/2.
qy_bin_centers = (qy_bin[1:]+qy_bin[:-1])/2.

dqx_bin_centers = (dqx_bin[1:]+dqx_bin[:-1])/2.
dqy_bin_centers = (dqy_bin[1:]+dqy_bin[:-1])/2.

# Grids for I, dqx, dqy
qx_bin_grid, qy_bin_centers_grid = np.meshgrid(qx_bin, qy_bin_centers)

i = np.sqrt((qx_bin_grid-0.45)**2+(qy_bin_centers_grid-0.35)**2) * 100
i_sigma = np.sqrt(i)


dqx_bin_centers_grid = np.tile(dqx_bin, (len(dqy_bin_centers), 1))
dqy_bin_centers_grid = np.tile(np.array([dqy_bin_centers]).transpose(),
                               (1, len(dqx_bin)))

# Q WS
iqxqy_ws2d = CreateWorkspace(
    DataX=qx_bin_grid,  # 2D
    DataY=i,  # 2D
    DataE=i_sigma,  # 2D
    NSpec=256,
    UnitX='MomentumTransfer',
    VerticalAxisUnit='MomentumTransfer',
    VerticalAxisValues=qy_bin_centers,  # 1D
)

dqx_ws2d = CreateWorkspace(
    DataX=qx_bin_grid,  # 2D
    DataY=dqx_bin_centers_grid,  # 2D
    DataE=None,  # 2D
    NSpec=256,
    UnitX='MomentumTransfer',
    VerticalAxisUnit='MomentumTransfer',
    VerticalAxisValues=qy_bin_centers,  # 1D
)

dqy_ws2d = CreateWorkspace(
    DataX=qx_bin_grid,  # 2D
    DataY=dqy_bin_centers_grid,  # 2D
    DataE=None,  # 2D
    NSpec=256,
    UnitX='MomentumTransfer',
    VerticalAxisUnit='MomentumTransfer',
    VerticalAxisValues=qy_bin_centers,  # 1D
)

GroupWorkspaces(InputWorkspaces=[iqxqy_ws2d, dqx_ws2d, dqy_ws2d],
                OutputWorkspace="qxqy_ws2d")

# Rebin the grouped WS
Rebin2D(
    InputWorkspace='qxqy_ws2d', OutputWorkspace='tmp',
    Axis1Binning='0.3,0.03,0.8', Axis2Binning='0.4,0.03,0.7')