====== labelThresholdEvents.m ======
**Class:** Function
**Description:** Takes in the [[usage:tracking#outputs|procTracks]] data structure and detects events based on the times at which the values of a user-specified field rise above (or drops below) a user-specified threshold. Applies a refractory period, preventing new events from being detected before the given quantity drops below (or rises above) a second user-specified threshold. This can help to suppress spurious events that result from noisy data.
**Author:** Oliver J. Meacock
function procTracks = labelThresholdEvents(procTracks,tgtPop,tgtField,startThresh,endThresh,eventLabel)
%LABELTHRESHOLDEVENTS labels times at which a given track-associated
%quantity either goes above or below a given threshold with an integer
%label.
%
% INPUTS:
% -procTracks: Track data output by the tracking module.
% -tgtPop: If population assignment has been performed, this variable
% can be used to specify which of the populations you wish to
% restrict event detection to. Can be left empty if you wish to apply
% it to all tracks.
% -tgtField: String defining the name of the field of procTracks you
% wish to apply threshold detection to. Can refer to any field,
% provided it is indexed by time and is one-dimensional (so 'x' and
% 'vmag' are permissible, but 'start' and 'length' are not).
% -startThresh: The value of the track-associated data that needs to
% be breached for an event to be detected. Track then enters a
% 'refractory' state, where the value needs to drop above/below
% endThresh for a new event to be detected.
% -endThresh: Defines the value at which the track-associated
% quantity is released from its refractory state. If endThresh >
% startThresh, it is assumed that an event starts when the given data
% drops below startThresh and that the refractory period ends when it
% rises above endThresh. If endThresh < startThresh, the opposite is
% true.
% -eventLabel: The integer you wish to use to label this set of
% events.
%
% OUTPUTS:
% -procTracks: Equal to the input procTracks, but either with the
% 'event' field updated, or (if not already present) added, with the
% appropriate events marked.
%
% Author: Oliver J. Meacock, (c) 2020
if startThresh >= endThresh
transType = 'up'; %Events are when you go from low to high
elseif startThresh < endThresh
transType = 'down'; %Events are when you go from high to low
else
error('Thresholds do not appear to be numerical.')
end
%Needed to ensure you don't overwrite pre-existing event data in procTracks
%later
if isfield(procTracks,'event')
preFlag = true;
else
preFlag = false;
end
for i = 1:size(procTracks,2)
currList = procTracks(i).(tgtField);
if preFlag
currEvt = procTracks(i).event; %Get pre-existing event list for this track (if it exists)
else
currEvt = zeros(size(procTracks(i).x)); %Want to make sure we're using the maximal possible length of the event list (need to avoid using e.g. vmag, which is one shorter than maximum as it is from a differential).
end
%Want to do processing only if the tgtPop variable is empty OR there is
%a target population, the population field exists, and the current
%track is part of the target population
if isempty(tgtPop)
proceed = true;
elseif ~isempty(tgtPop) && isfield(procTracks,'population')
if procTracks(i).population == tgtPop
proceed = true;
else
proceed = false;
end
elseif ~isempty(tgtPop) && ~isfield(procTracks,'population')
warning('Population specified, but field not present in procTracks. Proceeding without splitting by population...')
proceed = true;
end
if proceed
currInd = 1;
refState = false; %Whether you are in the refractory period or not
switch transType
case 'up'
evtSwitch = find(diff(currList > startThresh) == 1);
refSwitch = find(diff(currList < endThresh) == 1);
case 'down'
evtSwitch = find(diff(currList < startThresh) == 1);
refSwitch = find(diff(currList > endThresh) == 1);
end
evtList = [];
while proceed %I.e. while there are still events to test
futureEvts = [0,evtSwitch' > currInd];
futureRefs = [0,refSwitch' > currInd];
if ~refState && sum(futureEvts) > 0 %If event is permitted, detect next and switch to refractory state
currInd = evtSwitch(diff(futureEvts) == 1);
evtList = [evtList;currInd];
refState = true;
if sum(refSwitch' > currInd) == 0
proceed = false;
end
elseif refState && sum(futureRefs) > 0 %If event is not permitted, detect next refractory period end and switch to 'excitable' state
currInd = refSwitch(diff(futureRefs) == 1);
refState = false;
if sum(evtSwitch' > currInd) == 0
proceed = false;
end
else
proceed = false;
end
end
currEvt(evtList+1) = eventLabel;
end
procTracks(i).event = currEvt;
end
**Example usage:**
%Make some artificial data with multiple peaks of different heights
dt = 0.01;
t = 0:dt:20;
y = sin(x).*sin(x/4);
procTracks(1).times = round(t/dt);
procTracks(1).y = y;
%Detect event class 1 when y rises above -0.2 (must drop below -0.6 for
%subsequent detections)
procTracks = labelThresholdEvents(procTracks,[],'y',-0.2,-0.6,1);
%Detect event class 2 when y drops below -0.2 (must rise above 0.05 for
%subsequent detections
procTracks = labelThresholdEvents(procTracks,[],'y',-0.2,0.05,2);
%Plot results
plot(procTracks.times*dt,procTracks.y,'LineWidth',1.5)
hold on
%Red dots indicate events of class 1
plot(procTracks.times(procTracks(1).event == 1)*dt,y(procTracks(1).event == 1),'r.','MarkerSize',12)
%Black dots indicate events of class 2
plot(procTracks.times(procTracks(1).event == 2)*dt,y(procTracks(1).event == 2),'k.','MarkerSize',12)
**Example output:**
{{ :post:eventlabelexample.png?nolink&400 |}}