Merge branch 'devel' of iffgit.fz-juelich.de:zolnie/pysen into devel

def plot_curscan(axes, filename, ring=None,

                 normalize=False, total=False, legend=True):

def plot_curscan(axes, filename, ring=None, normalize=False,

                 currents=False, total=False, legend=True):

    "plot current scan"

    hdr, res = read_datfile(filename)

    #q     = hdr['q'][0]

              hdr['phase_points_spindown'] + hdr['phase_points_spinup'] )

    phases = []

    _phase = []

    cur   = []

    for l in res:

        #print "========================================="

        pc = 0

                pc = v

            if k.startswith('detector.ring') and k.endswith('sum'):

                _phase.append([v, pc])

            if k == 'phasecurrent1':

                cur.append(v)

        phases.append(_phase)

    phases = np.asarray(phases)

    cur    = np.asarray(cur)

    def plot_counts(counts, label, fmt='ko-', normalize=None):

        "plot counts as a function of phase point"

        iphase  = range(1, len(counts)+1)

        if currents:

            x= cur

        else:

            x = range(1, len(counts)+1)

        ecounts = np.sqrt(counts*1.0)

        if normalize is not None:

            counts  = counts *1.0/normalize

            ecounts = ecounts*1.0/normalize

        axes.errorbar(iphase, counts, yerr=ecounts,

        axes.errorbar(x, counts, yerr=ecounts,

                      fmt=fmt, label=label) # 'Ring #%d' % (r+1))

    if total:

    if legend:

        axes.legend(loc='upper right')

        if not currents:

            axes.set_xlim((0, ntot+8))

        axes.set_title(os.path.basename(hdr['info'].split()[-1]))

    else:

        if not currents:

            axes.set_xlim((0, ntot+2))

    axes.grid(True)

#!/usr/bin/env python

"ATARI Plot"

import numpy as np

import matplotlib.pyplot as plt

import h5py

from matplotlib import cm

from functools import partial

from pysen import ANGSTROM, NANOSECOND, OMEGA_N

from pysen.mathutil import normalize_counts

from pysen.echo.fit import flux_weighted_eshape, echo_fit, echo_curve

def atari_plot(hdfile, iecho=0, tbin1=0, tbin2=None, npix=1, iphase=15):

    "atari plot"

    ntaus   =  hdfile.attrs['scan_length']

    n_idx = { 'dn': hdfile['/'].attrs['point_to_down'],

              'up': hdfile['/'].attrs['point_to_up'] ,

              'nphases': hdfile['/'].attrs['no_of_phases'] }

    nt      = hdfile['/detector'].attrs['no_t_channels']

    ny      = hdfile['/detector'].attrs['no_y_channels']

    nx      = hdfile['/detector'].attrs['no_x_channels']

    nph     = n_idx['dn']

    for echo in list(hdfile['/data'].values()):

        if iecho is not None and echo.attrs['id'] != iecho:

            print("skipping", echo.name)

            continue

        else:

            print("processing", echo.name)

        phase   = echo['phase']

        physics = echo['phys']

        params  = echo['params']

        tau0    = physics.attrs['fouriertime']/NANOSECOND

        j0      = physics.attrs['prim_fieldintegral']

        lam0    = physics.attrs['lambda']

        phase0  = float(echo['tech'].attrs['i5'][0])

        lmax    = np.array(params['lambdaTable']).flatten()/ANGSTROM

        dlam    = np.ones_like(lmax)*(lmax[1]-lmax[0])

        tau     = np.array(params['fouriertimeTable']).flatten()/NANOSECOND

        tau     = np.tile(tau, (ny, nx)).reshape((ny, nx, nt))

        det  = phase['detector'][...]

        moni = phase['counter' ][...]

        pcha = phase['proton_charge'][...]

        cur  = phase['phase_current'][:, 0] # actual value

        pcha0 = np.average(pcha)

        pcha = np.tile(pcha, (ny, nx, nt)).reshape(n_idx['nphases'],ny,nx,nt)

        pcha = pcha/float(pcha0)

        det  = det/pcha

        #

        plt.subplot(2,2,1)

        wlen = np.sum(det, axis=(0,1,2))

        #plt.step(lmax[tbin1:tbin2], wlen[tbin1:tbin2])

        plt.step(wlen[tbin1:tbin2], lmax[tbin1:tbin2])

        plt.ylabel(r'$\lambda$ [$\AA$]')

        plt.xlim(xmin=0)

        plt.grid(True)

        plt.subplot(2,2,2)

        xatari = np.sum(det, axis=(1,2)).T

        extent=(cur[0], cur[nph-1], lmax[tbin1], lmax[tbin2])

        print(cur[nph])

        plt.imshow(xatari[tbin1:tbin2,:nph], origin='lower', cmap=cm.hot, extent=extent, aspect='auto')

        plt.subplot(2,2,3)

        img = np.sum(det[:,:,:,tbin1:tbin2], axis=(0,-1))

        plt.imshow(img , cmap=cm.jet, aspect='auto')

        plt.subplot(2,2,4)

        pha = np.sum(det[:nph,:,:,tbin1:tbin2], axis=(1,2,3))

        eph = np.sqrt(pha)

        imx = iphase or np.argmax(pha)

        deg = 4

        dix = 4

        imx1 = max(0      ,  imx-dix)

        imx2 = min(len(pha), imx+dix)

        xc = cur[imx1:imx2+1]

        yp = pha[imx1:imx2+1]

        pf = np.polyfit(xc, yp, deg=deg)

        xf = np.linspace(min(xc), max(xc), 1001)

        yf = np.polyval(pf, xf)

        imx = np.argmax(yf)

        x_fit = xf[imx]

        plt.errorbar(cur[:nph], pha, yerr=eph, fmt='ks--')

        plt.plot(xf, np.polyval(pf, xf), 'r-') #, label='fit (poly%d)' % deg)

        plt.axvline(x_fit, color='green', label='fit=%.3f A' % x_fit)

        plt.axvline(phase0, color='blue')

        plt.xlabel(r'phase current [A]')

        plt.legend(loc='best')

        plt.grid(True)

def main():

    "the main"

    import argparse

    parser = argparse.ArgumentParser(description='atari plot')

    parser.add_argument('file', metavar='filename', help='file to process', nargs='+')

    parser.add_argument('--echo',  '-e',  dest='echo', type=int,   default=None,

                        help='set echo to display (default=%(default)s)')

    parser.add_argument('--t1',  '-t',  dest='tbin1', type=int,   default=1,

                        help='set min TOF bin (default=%(default)s)')

    parser.add_argument('--t2',  '-T',  dest='tbin2', type=int,   default=-1,

                        help='set max TOF bin (default=%(default)s)')

    parser.add_argument('--phase','-p',  dest='phase', type=int,   default=None,

                        help='set center phase(default=%(default)s)')

    #parser.add_argument('--npix', '-n', dest='n', type=int, default=1,

    #                    help='sum NxN pixels')

    args = parser.parse_args()

    for filename in args.file:

        with h5py.File(filename, 'r') as hdf5file:

            plt.figure(figsize=(8,8))

            atari_plot(hdf5file, iecho=args.echo, tbin1=args.tbin1, tbin2=args.tbin2, iphase=args.phase)

    plt.show()

if __name__ == "__main__":

    main()

    parser = argparse.ArgumentParser(description=description)

    parser.set_defaults(normalize=False, total=False, out=None, ring=None,

                        legend=True, show=True, live=None,

                        legend=True, show=True, live=None, currents=False,

                        sx=8.0, sy=8.0, dpi=100)

    parser.add_argument('file', metavar='inpfile',

                        help='file to process')

                        help='plot total counts')

    parser.add_argument('--ring', '-r', dest='ring', type=int,

                        help='plot total counts')

    parser.add_argument('--current', '-C', dest='currents', action='store_true',

                        help='plot current instead of point index')

    parser.add_argument('--out', '-o', dest='out',

                        help='set outfile title')

    parser.add_argument('--live', '-i', dest='live', metavar='sec', type=int,

        print(i, args.file)

        qpl.plot_curscan(axes, args.file,

                         normalize=args.normalize, ring=args.ring,

                         total=args.total, legend=args.legend)

                         total=args.total, currents=args.currents, legend=args.legend)

    if args.live:

        _ = animation.FuncAnimation(fig, plot_curscan, interval=int(args.live)*1000)

        'scripts/filemon.py',

        'scripts/qbin.py',

        'scripts/qtau.py',

        'scripts/atariplot.py',

        #'scripts/echoplot.py',

        #'scripts/stapler.py',

        #'scripts/the_cube.py',