Saccade detection

import numpy as np
from matplotlib import pyplot as plt

# Source files for data
src = '/home/sebastiaan/Data/0045 Pupil cuing/0045E Pre-saccadic pupil/traces/0045E00A.asc-%.5d-%s.npy'

vThr = 2 # Velocity threshold for both saccade onset and offset
aThr = .1 # Acceleration threshold for saccade onset
wLen = 5 # Moving window size. Larger values are more accurate but slower.

for trialId in range(100):

    # Load the data into a numpy array with three columns: x, y, pupil size
    a = np.concatenate([
        np.load(src % (trialId, 'baseline')),
        np.load(src % (trialId, 'cue')),
        np.load(src % (trialId, 'postSacc'))
        ])


    # We will use this array to maintain the eye position during the moving
    # window average.
    aPos = np.zeros( [wLen, 2] )

    # These arrays are to store the outcome so that we can plot it later on.
    # They don't play a role in the detection algorithm.
    _vel = np.zeros(len(a))
    _acc = np.zeros(len(a))
    _sacc = np.zeros(len(a))

    inSacc = 0 # 0 or 1 depending on whether we are in a saccade

    for i in range(len(a)):

        # Get a sample. Here we read the x, y coordinates from the data.
        # Generally you would query the eye tracker.
        x = a[i,0]
        y = a[i,1]

        # Shift the array one step forward and insert the new sample as the
        # first entry in the array.
        aPos[1:] = aPos[:-1]
        aPos[0] = [x,y]

        # Get the coordinates difference scores for consecutive samples. This is
        # a 2D array (so velocity for x and y separately).
        d = aPos[1:] - aPos[:-1]
        # Calculate the velocity. This is a 1D array, which is the Pythagoras of
        # the previous array.
        aVel = np.sqrt(d[:,0]**2 + d[:,1]**2)
        # Calculate the difference scores for the velocity, i.e. the
        # acceleration.
        aAcc = aVel[:-1] - aVel[1:]

        # Calculate the mean velocity and acceleration over the moving window.
        mVel = aVel.mean()
        mAcc = aAcc.mean()

        # If either the velocity or the acceleration threshold is exceeded,
        # detect a saccade.
        if mVel > vThr or mAcc > aThr:
            inSacc = 1
        # If the velocity drops back under threshold (you can use a lower
        # threshold here as well, may be) detect and-end of-saccade.
        if mVel < vThr:
            inSacc = 0

        # Save for plotting
        _vel[i] = mVel
        _acc[i] = mAcc
        _sacc[i] = inSacc

    # Plot to see it worked!
    yPos = a[:,0]
    plt.plot(yPos[wLen:] / 100, color='green')
    plt.plot(_vel[wLen:], color='red')
    plt.plot(_acc[wLen:], color='blue')
    plt.plot(_sacc[wLen:] * 10, color='black')
    plt.show()