11"""
2- Set of functions to handle wheel data
2+ Set of functions to handle wheel data.
33"""
4+ import logging
5+ import warnings
6+ import traceback
7+
48import numpy as np
59from numpy import pi
10+ from iblutil .numerical import between_sorted
611import scipy .interpolate as interpolate
712import scipy .signal
813from scipy .linalg import hankel
1318__all__ = ['cm_to_deg' ,
1419 'cm_to_rad' ,
1520 'interpolate_position' ,
16- 'last_movement_onset ' ,
21+ 'get_movement_onset ' ,
1722 'movements' ,
1823 'samples_to_cm' ,
1924 'traces_by_trial' ,
20- 'velocity_smoothed ' ]
25+ 'velocity_filtered ' ]
2126
2227# Define some constants
2328ENC_RES = 1024 * 4 # Rotary encoder resolution, assumes X4 encoding
@@ -49,6 +54,8 @@ def interpolate_position(re_ts, re_pos, freq=1000, kind='linear', fill_gaps=None
4954 Timestamps of interpolated positions
5055 """
5156 t = np .arange (re_ts [0 ], re_ts [- 1 ], 1 / freq ) # Evenly resample at frequency
57+ if t [- 1 ] > re_ts [- 1 ]:
58+ t = t [:- 1 ] # Occasionally due to precision errors the last sample may be outside of range.
5259 yinterp = interpolate .interp1d (re_ts , re_pos , kind = kind )(t )
5360
5461 if fill_gaps :
@@ -63,7 +70,7 @@ def interpolate_position(re_ts, re_pos, freq=1000, kind='linear', fill_gaps=None
6370
6471def velocity (re_ts , re_pos ):
6572 """
66- Compute wheel velocity from non-uniformly sampled wheel data. Returns the velocity
73+ (DEPRECATED) Compute wheel velocity from non-uniformly sampled wheel data. Returns the velocity
6774 at the same samples locations as the position through interpolation.
6875
6976 Parameters
@@ -78,6 +85,13 @@ def velocity(re_ts, re_pos):
7885 np.ndarray
7986 numpy array of velocities
8087 """
88+ for line in traceback .format_stack ():
89+ print (line .strip ())
90+
91+ msg = 'brainbox.behavior.wheel.velocity has been deprecated. Use velocity_filtered instead.'
92+ warnings .warn (msg , DeprecationWarning )
93+ logging .getLogger (__name__ ).warning (msg )
94+
8195 dp = np .diff (re_pos )
8296 dt = np .diff (re_ts )
8397 # Compute raw velocity
@@ -92,12 +106,23 @@ def velocity(re_ts, re_pos):
92106
93107def velocity_filtered (pos , fs , corner_frequency = 20 , order = 8 ):
94108 """
95- Compute wheel velocity from uniformly sampled wheel data
109+ Compute wheel velocity from uniformly sampled wheel data.
110+
111+ pos: array_like
112+ Vector of uniformly sampled wheel positions.
113+ fs : float
114+ Frequency in Hz of the sampling frequency.
115+ corner_frequency : float
116+ Corner frequency of low-pass filter.
117+ order : int
118+ Order of Butterworth filter.
96119
97- :param pos: vector of uniformly sampled wheel positions
98- :param fs: scalar, sampling frequency
99- :param corner_frequency: scalar, corner frequency of low-pass filter
100- :param order: scalar, order of Butterworth filter
120+ Returns
121+ -------
122+ vel : np.ndarray
123+ Array of velocity values.
124+ acc : np.ndarray
125+ Array of acceleration values.
101126 """
102127 sos = scipy .signal .butter (** {'N' : order , 'Wn' : corner_frequency / fs * 2 , 'btype' : 'lowpass' }, output = 'sos' )
103128 vel = np .insert (np .diff (scipy .signal .sosfiltfilt (sos , pos )), 0 , 0 ) * fs
@@ -107,7 +132,7 @@ def velocity_filtered(pos, fs, corner_frequency=20, order=8):
107132
108133def velocity_smoothed (pos , freq , smooth_size = 0.03 ):
109134 """
110- Compute wheel velocity from uniformly sampled wheel data
135+ (DEPRECATED) Compute wheel velocity from uniformly sampled wheel data.
111136
112137 Parameters
113138 ----------
@@ -125,6 +150,13 @@ def velocity_smoothed(pos, freq, smooth_size=0.03):
125150 acc : np.ndarray
126151 Array of acceleration values
127152 """
153+ for line in traceback .format_stack ():
154+ print (line .strip ())
155+
156+ msg = 'brainbox.behavior.wheel.velocity_smoothed has been deprecated. Use velocity_filtered instead.'
157+ warnings .warn (msg , DeprecationWarning )
158+ logging .getLogger (__name__ ).warning (msg )
159+
128160 # Define our smoothing window with an area of 1 so the units won't be changed
129161 std_samps = np .round (smooth_size * freq ) # Standard deviation relative to sampling frequency
130162 N = std_samps * 6 # Number of points in the Gaussian covering +/-3 standard deviations
@@ -141,15 +173,24 @@ def velocity_smoothed(pos, freq, smooth_size=0.03):
141173
142174def last_movement_onset (t , vel , event_time ):
143175 """
144- Find the time at which movement started, given an event timestamp that occurred during the
145- movement. Movement start is defined as the first sample after the velocity has been zero
146- for at least 50ms. Wheel inputs should be evenly sampled.
176+ (DEPRECATED) Find the time at which movement started, given an event timestamp that occurred during the
177+ movement.
178+
179+ Movement start is defined as the first sample after the velocity has been zero for at least 50ms.
180+ Wheel inputs should be evenly sampled.
147181
148182 :param t: numpy array of wheel timestamps in seconds
149183 :param vel: numpy array of wheel velocities
150184 :param event_time: timestamp anywhere during movement of interest, e.g. peak velocity
151185 :return: timestamp of movement onset
152186 """
187+ for line in traceback .format_stack ():
188+ print (line .strip ())
189+
190+ msg = 'brainbox.behavior.wheel.last_movement_onset has been deprecated. Use get_movement_onset instead.'
191+ warnings .warn (msg , DeprecationWarning )
192+ logging .getLogger (__name__ ).warning (msg )
193+
153194 # Look back from timestamp
154195 threshold = 50e-3
155196 mask = t < event_time
@@ -166,6 +207,42 @@ def last_movement_onset(t, vel, event_time):
166207 return t
167208
168209
210+ def get_movement_onset (intervals , event_times ):
211+ """
212+ Find the time at which movement started, given an event timestamp that occurred during the
213+ movement.
214+
215+ Parameters
216+ ----------
217+ intervals : numpy.array
218+ The wheel movement intervals.
219+ event_times : numpy.array
220+ Sorted event timestamps anywhere during movement of interest, e.g. peak velocity, feedback
221+ time.
222+
223+ Returns
224+ -------
225+ numpy.array
226+ An array the length of event_time of intervals.
227+
228+ Examples
229+ --------
230+ Find the last movement onset before each trial response time
231+
232+ >>> trials = one.load_object(eid, 'trials')
233+ >>> wheelMoves = one.load_object(eid, 'wheelMoves')
234+ >>> onsets = last_movement_onset(wheelMoves.intervals, trials.response_times)
235+ """
236+ if not np .all (np .diff (event_times ) > 0 ):
237+ raise ValueError ('event_times must be in ascending order.' )
238+ onsets = np .full (event_times .size , np .nan )
239+ for i in np .arange (intervals .shape [0 ]):
240+ onset = between_sorted (event_times , intervals [i , :])
241+ if np .any (onset ):
242+ onsets [onset ] = intervals [i , 0 ]
243+ return onsets
244+
245+
169246def movements (t , pos , freq = 1000 , pos_thresh = 8 , t_thresh = .2 , min_gap = .1 , pos_thresh_onset = 1.5 ,
170247 min_dur = .05 , make_plots = False ):
171248 """
@@ -296,7 +373,7 @@ def movements(t, pos, freq=1000, pos_thresh=8, t_thresh=.2, min_gap=.1, pos_thre
296373 if make_plots :
297374 fig , axes = plt .subplots (nrows = 2 , sharex = 'all' )
298375 indices = np .sort (np .hstack ((onset_samps , offset_samps ))) # Points to split trace
299- vel , acc = velocity_smoothed (pos , freq , 0.015 )
376+ vel , acc = velocity_filtered (pos , freq )
300377
301378 # Plot the wheel position and velocity
302379 for ax , y in zip (axes , (pos , vel )):
@@ -440,6 +517,6 @@ def to_mask(a, b):
440517 return [(cuts [n ][0 , :], cuts [n ][1 , :]) for n in range (len (cuts ))] if separate else cuts
441518
442519
443- if __name__ == " __main__" :
520+ if __name__ == ' __main__' :
444521 import doctest
445522 doctest .testmod ()
0 commit comments