function [tm, iPre, iPost, iLast] ...
  = getTransitions(x, t, upperBound, lowerBound, upperRef, midRef, lowerRef)
%GETTRANSITIONS Extract indices and instants that bound transitions.
% This function finds the transitions of the waveform between the specified
% lower and upper reference levels as described in IEEE Std 181-2003 Sec
% 5.3 and returns results with respect to positive and negative going
% transitions.
% 
% For each identified transition it returns the transition instants
% that correspond to the 50% crossing between level1 and level2 and the
% corresponding transition durations (either rising or falling) where the
% signal crosses the x1 and x2 parameters.
%
% Input Arguments:
%  x           a column vector containing the sample amplitude values
%  t           a vector containing the corresponding sample instants.
%
%  upperBound  the state boundaries of the upper and lower level
%  lowerBound  used to detect transitions
%
%  lowerRef    The reference levels that bound a transition.   
%  upperRef    (These are typically the 10% and 90% reference levels).
%
%  midRef      The reference level used to determine the transition instant.
%              This is typically the 50% reference level.
%
% Output Arguments:
%
%  tm       A structure with the following fields of column vectors:
%
%      Polarity    the polarities of each transition
%                  (+1) for a rising  (positive going)
%                  (-1) for a falling (negative-going) transition.
%
%      Duration    The transition durations that correspond to each
%                  transition.  These are computed from the time it takes
%                  to cross the lower and upper reference level.
%
%      LowerCross  The instants where each transitions crosses the lower
%                  reference level
%
%      MiddleCross The transition instants that correspond to where the
%                  signal crosses the 50% point between level1 and level2
%
%      UpperCross  The instant where each transition crosses the upper
%                  reference level
%
%      SlewRate    the ratio of (upperCrossing - lowerCrossing)/duration
%
%  iPre     The indices just before each transition begins.
%
%  iPost    The indices just after each transition ends.
%  
%  iLast    The index of the last sample known to be in either state.

%   Copyright 2011-2013 The MathWorks, Inc.
%#codegen
   
  %find indices located outside a transition region
  iState     = find(x<lowerBound | x>upperBound);
  if numel(iState) > 1
    iStatePre  = iState(1:(end-1));
    iStatePost = iState(2:end);
        
    % look for positive going and negative going transitions.
    iTrans  = find(x(iStatePre) < lowerBound & upperBound < x(iStatePost) | ...
                   x(iStatePre) > upperBound & lowerBound > x(iStatePost) );
  
    iPre    = iStatePre(iTrans);
    iPost   = iStatePost(iTrans);

    % Mark the polarities of each transition
    polarity = 2 * (x(iPre) < lowerBound) - 1;

    % Set the reference levels for the pre and post transitions
    preRef  = signal.internal.polaritySelect(polarity, lowerRef, upperRef);
    postRef = signal.internal.polaritySelect(polarity, upperRef, lowerRef);

    %Now loop through and extract the transition occurrence instants
    %see IEEE Standard 181-2003 section 5.3.3.1
    numTrans = numel(iTrans);
    iRMid  = zeros(numTrans, 1);
    iRPre  = zeros(numTrans, 1);
    iRPost = zeros(numTrans, 1);
  else
    iPre = zeros(0,1);
    iPost = zeros(0,1);
    polarity = zeros(0,1);
    preRef = zeros(0,1);
    postRef = zeros(0,1);
    numTrans = 0;
    iRMid = zeros(0,1);
    iRPre = zeros(0,1);
    iRPost = zeros(0,1);
  end
  
  for i = 1:numTrans
    %define convenience indices for compactness 
    iA = iPre(i);
    iB = iPost(i);

    if polarity(i) > 0
      %look for the first positive crossing of midRef (5.3.3.1)
      iX = find((x(iA:iB-1) <= midRef & midRef < x(iA+1:iB)), 1, 'first');
      iRMid(i) = iX + iA - 1;

      % look for the positive preRef crossing closest to iRMid(i)
      iX = find(x(iA:iRMid(i)+1) < preRef(i), 1, 'last');
      iRPre(i) = iX + iA - 1;

      % look for the positive postRef crossing closest to iRMid(i)
      iX = find(postRef(i) < x(iRMid(i):iB), 1, 'first');
      iRPost(i) = iX + iRMid(i) - 2;
    else
      %look for the first negative crossing of midRef (5.3.3.1)
      iX = find(x(iA:iB-1)>= midRef & midRef > x(iA+1:iB), 1, 'first');
      iRMid(i) = iX + iA - 1;

      % look for the negative preRef crossing closest to iRMid(i)
      iX = find(x(iA:iRMid(i)+1) >= preRef(i), 1, 'last');
      iRPre(i) = iX + iA - 1;

      % look for the negative postRef crossing closest to iRMid(i)
      iX = find(postRef(i) > x(iRMid(i):iB), 1, 'first');
      iRPost(i) = iX + iRMid(i) - 2;
    end
    
  end

  if numTrans > 0
    %interpolate to find the corresponding reference level instant
    tMiddle = signal.internal.linterp(t(iRMid), t(iRMid+1), x(iRMid), x(iRMid+1), midRef);
    tRPre   = signal.internal.linterp(t(iRPre), t(iRPre+1), x(iRPre), x(iRPre+1), preRef);
    tRPost  = signal.internal.linterp(t(iRPost),t(iRPost+1),x(iRPost),x(iRPost+1),postRef);

    %Now extract the corresponding durations
    tDuration = tRPost - tRPre;

    %Now compute the slew rate
    rSlew = (postRef - preRef) ./ tDuration;

    %Now obtain the upper and lower crossing instants
    tUpper = signal.internal.polaritySelect(polarity, tRPost, tRPre);
    tLower = signal.internal.polaritySelect(polarity, tRPre, tRPost);
  else
    tMiddle   = zeros(0,1);
    tDuration = zeros(0,1);
    rSlew     = zeros(0,1);
    tUpper    = zeros(0,1);
    tLower    = zeros(0,1);
    polarity  = zeros(0,1);
  end

  
  %Now extract the index of the last known state.  This will help us know
  %how much of the waveform needs to be saved when seeking the next
  %transition.
  if numel(iState) > 0
    iLast = iState(end);
  else
    iLast = 0;  %no state found
  end
  
  % Now build the transition metrics structure
  tm = makeTransitionMetrics(tDuration, polarity, rSlew, tMiddle, ...
                             tLower, tUpper);
end

% transition metrics
function tm = ...
  makeTransitionMetrics(tDuration, polarity, rSlew, tMiddle, tLower, tUpper)
  tm = struct('Duration',    tDuration, ...
              'Polarity',    polarity, ...
              'SlewRate',    rSlew, ...
              'MiddleCross', tMiddle, ...
              'LowerCross',  tLower, ...
              'UpperCross',  tUpper);
end
