Commit b022b990 authored by Zolnierczuk, Piotr's avatar Zolnierczuk, Piotr
Browse files

updates to xyz plot

parent 11f552eb
Loading
Loading
Loading
Loading
+79 −59
Original line number Diff line number Diff line
@@ -20,6 +20,52 @@ def _norm_pix(x1, x2, nx, whole_detector=False):
    return np.mod(x1,nx), np.mod(x2,nx)


def xyz_decomposition(xyzdata):
    """Decompose {x,y,z}_{up,down} measurements
    into nuclear coherent, incoherent and magnetic components

    Data is expected to be a dictionary of numbers or numpy arrays
    with six keys: {x,y,z}_{up,dn}"""

    res = {}
    # sum of x, y and z  measurements
    sX  = xyzdata['x_up'] + xyzdata['x_dn']
    sY  = xyzdata['y_up'] + xyzdata['y_dn']
    sZ  = xyzdata['z_up'] + xyzdata['z_dn']

    # diff up-down for x, y and z measurements
    dX  = xyzdata['x_up'] - xyzdata['x_dn']
    dY  = xyzdata['y_up'] - xyzdata['y_dn']
    dZ  = xyzdata['z_up'] - xyzdata['z_dn']

    # sum of all up
    sUP    =    xyzdata['z_up'] + xyzdata['x_up'] + xyzdata['y_up']
    sDN    =    xyzdata['z_dn'] + xyzdata['x_dn'] + xyzdata['y_dn']
    # "magnetic" sum of up/down
    sM_UP  = +2*xyzdata['z_up'] - xyzdata['x_up'] - xyzdata['y_up']
    sM_DN  = -2*xyzdata['z_dn'] + xyzdata['x_dn'] + xyzdata['y_dn']

    # the results
    sTOT= (sUP + sDN)/3     # total "cross section"
    sM  =  2*dZ - (dX + dY) # magnetic "cross section"
    sN  = (2*sUP - sDN)/6   # nuclear coherent
    sI  = sTOT - sN - sM    # incoherent

    #
    res['n_coh'] = sN
    res['i_inc'] = sI
    res['m_mag'] = sM
    #
    res['sum_ave'] = sTOT
    res['sum_x'  ] = sX
    res['sum_y'  ] = sY
    res['sum_z'  ] = sZ
    res['m_ave'  ] = sM_UP + sM_DN
    res['m_up/2']  = sM_UP
    res['m_dn/2']  = sM_DN
    return res


def plot_xyz_single(hdfile, fig, axis, axis1=None, **kwargs):
    "xyz plot"
    #
@@ -139,7 +185,7 @@ def plot_xyz_single(hdfile, fig, axis, axis1=None, **kwargs):



def xyz_analysis(results, axis, norm=False, output=None, details=False):
def xyz_analysis(results, fig, axes, norm=False, output=None, details=False):
    "XYZ decomposition"

    pcha = []
@@ -150,53 +196,38 @@ def xyz_analysis(results, axis, norm=False, output=None, details=False):
        pcha = 1
    cnorm    = {}
    cdetails = {}
    ax_xyz, axn, axi, axm = axes  # xyz, nucl, inc, mag

    if fig is None:
        fig = plt.gcf()

    for key, val in results.items():
        cnorm[key] = val['counts']/val['proton_charge']*pcha

    sX  = cnorm['x_up'] + cnorm['x_dn']
    sY  = cnorm['y_up'] + cnorm['y_dn']
    sZ  = cnorm['z_up'] + cnorm['z_dn']
    sUP = cnorm['x_up'] + cnorm['y_up'] + cnorm['z_up']
    sDN = cnorm['x_dn'] + cnorm['y_dn'] + cnorm['z_dn']
    #
    dX  = cnorm['x_up'] - cnorm['x_dn']
    dY  = cnorm['y_up'] - cnorm['y_dn']
    dZ  = cnorm['z_up'] - cnorm['z_dn']
    #
    sM_UP  = +2*cnorm['z_up'] - (cnorm['x_up'] + cnorm['y_up'])
    sM_DN  = -2*cnorm['z_dn'] + (cnorm['x_dn'] + cnorm['y_dn'])
    sTOT= (sUP + sDN)/3
    sM  =  2*dZ - (dX + dY)
    sN  = (2*sUP - sDN)/6
    sI  = sTOT - sN - sM

    cnorm['m_mag'] = sM
    cnorm['n_coh'] = sN
    cnorm['i_inc'] = sI

    cdetails['sum_ave'] = sTOT
    cdetails['sum_x'  ] = sX
    cdetails['sum_y'  ] = sY
    cdetails['sum_z'  ] = sZ
    #cdetails['m_ave'  ] = sM_UP + sM_DN
    cdetails['m_up/2']  = sM_UP
    cdetails['m_dn/2']  = sM_DN
    xyz = xyz_decomposition(cnorm)

    for lbl in ('n_coh', 'i_inc', 'm_mag') :
        cnorm[lbl] = xyz[lbl]

    for lbl in ('sum_ave', 'sum_x', 'sum_y', 'sum_z', 'm_up/2', 'm_dn/2'):
        cdetails[lbl] = xyz[lbl]

    sTOT = xyz['sum_ave']

    keys = list(cnorm.keys())
    vals = np.asarray(list(cnorm.values()))
    axis.step(keys, vals, where='mid')
    axis.grid(True)

    ax_xyz.step(keys, vals, where='mid')
    ax_xyz.grid(True)
    # fixing yticks with matplotlib.ticker "FixedLocator"
    axis.xaxis.set_ticks(range(len(keys)))
    axis.set_xticklabels(keys, rotation=40 , ha='right')
    ax_xyz.xaxis.set_ticks(range(len(keys)))
    ax_xyz.set_xticklabels(keys, rotation=40 , ha='right')
    ax_xyz.set_ylim(bottom=min(min(vals),0))

    out = StringIO()
    out.write("#label   counts \taverage q \taverage wavelength\n")
    for key, cval in cnorm.items():
        if key.startswith('m_mag'):
        if key.startswith('n_coh'):
            out.write("#summary\n")
        if norm:
            out.write("%-8.8s %10.6f" % (key, cval))
@@ -227,7 +258,14 @@ def xyz_analysis(results, axis, norm=False, output=None, details=False):
    else:
        print(out.getvalue())

    # summary plots
    image = { key : results[key]['image'] for key in results }
    imres = xyz_decomposition(image)

    for ax,key in zip((axn,axi,axm),('n_coh','i_inc','m_mag')):
        im = ax.imshow(imres[key], aspect='equal')
        fig.colorbar(im, ax=ax, use_gridspec=True)
        ax.set_title(key)

def plot_xyz(*hdfiles, **kwargs):
    "plot xyz"
@@ -255,6 +293,11 @@ def plot_xyz(*hdfiles, **kwargs):
    ax_lam = fig.add_subplot(gs1[1] )
    ax_q   = fig.add_subplot(gs1[2] )
    #
    gs2  = gs[2,:].subgridspec(1,3)
    ax_n = fig.add_subplot(gs2[0])
    ax_i = fig.add_subplot(gs2[1])
    ax_m = fig.add_subplot(gs2[2])
    #
    for hfilename in hdfiles:
        log.info('processing %s', hfilename)
        basename , _ = os.path.splitext(os.path.basename(hfilename))
@@ -273,32 +316,9 @@ def plot_xyz(*hdfiles, **kwargs):
                axis1=(ax_lam, ax_q)
            res = plot_xyz_single(hdf5file, fig, ax, axis1=axis1, label=label, **kwargs)
            results.update({label:res})
    #
    gsm = gs[2,0]
    axm = fig.add_subplot(gsm)

    dX  = results['x_up']['image'] + results['x_dn']['image']
    dY  = results['y_up']['image'] + results['y_dn']['image']
    dZ  = results['z_up']['image'] + results['z_dn']['image']
#    dX  = cnorm['x_up'] - cnorm['x_dn']
#    dY  = cnorm['y_up'] - cnorm['y_dn']
#    dZ  = cnorm['z_up'] - cnorm['z_dn']
#    #
#    sM_UP  = +2*cnorm['z_up'] - (cnorm['x_up'] + cnorm['y_up'])
#    sM_DN  = -2*cnorm['z_dn'] + (cnorm['x_dn'] + cnorm['y_dn'])
#    sTOT= (sUP + sDN)/3
    sM  =  2*dZ - (dX + dY)
    im = axm.imshow(sM , aspect='equal')
    fig.colorbar(im, ax=axm, use_gridspec=True)
#    sN  = (2*sUP - sDN)/6
#    sI  = sTOT - sN - sM

#    cnorm['m_mag'] = sM
#    cnorm['n_coh'] = sN
#    cnorm['i_inc'] = sI

    #
    xyz_analysis(results, ax_xyz, norm=normalize, output=savefile, details=details)
    xyz_analysis(results, fig, (ax_xyz,ax_n, ax_i, ax_m), norm=normalize, output=savefile, details=details)

    #
    title = title.replace('_', ' ')
+2 −2
Original line number Diff line number Diff line
@@ -3,8 +3,8 @@ PySEN revision module
"""
import sys
__version__  = "1.3"
__release__  = "a3"
__date__     = "June 1, 2023"
__release__  = "a4"
__date__     = "June 2, 2023"

def version(full=False):
    "get pysen version number"