Можно ли векторизовать эту функцию Matlab (или заставить работать быстрее другим методом)?

Question

Основная проблема заключается в поиске значений с фиксированным смещением от текущего значения.Мой текущий метод очень медленный, когда вектор Значения большой (обычно 100000 элементов).

function [ AverageValue ] = CalculateAverageValueOverAngle( Value, Angle )
% function [ AverageValue ] = CalculateAverageValueOverAngle( Value, Angle )
%   Calculate average value from instantaneous value and angle
%   Average value is calculated over +- 90 degrees from current angle

   AverageValue = zeros( size( Value ) );
   UnwrappedRadians = unwrap( Angle ./ 180 * pi );

   for i=1:length(UnwrappedRadians)
       mid = UnwrappedRadians(i);
       start = find( UnwrappedRadians(1:i) < (mid - pi/2), 1, 'Last');
       finish = find( UnwrappedRadians(i:end) > (mid + pi/2), 1, 'First');
       if isempty(start) | isempty(finish)
           AverageValue(i) = Value(i); 
        else
           AverageValue(i) = mean(Value(start:finish+i-1));  % nanmean
        end
    end    
end

Answer 1

MatlabSorter предоставил рефакторинг, который имеет смысл для меня, поэтому, если ваш код действительно делает то, что вы хотите, его рефакторинг - это путь вперед :-). Обратите внимание, что в моих тестах с numberl (Angle) = 50000 его рефакторинг не сильно сэкономил (возможно, потому что мои выборочные данные предполагали, что find () почти никогда не потерпит неудачу, за исключением самого начала и в конце трассировки данных) .

Однако, глядя на ваш код, мне стало интересно: вы уверены, что абсолютно хотите усреднить все значения с первого раза, когда угол попадает в диапазон mid-pi / 2 ... mid + pi / 2, до в последний раз он покидает этот диапазон? Если развернутые углы являются немонотонными (например, если разрешены обратные перемещения, если частота дискретизации слишком мала, чтобы избежать наложения, или просто из-за шума измерения в угле), вы также будете усреднять по некоторым значениям за пределами (и, возможно, за пределами) диапазона 180 °.

Обратите внимание, что в любом случае первое измеренное вами значение (значение (начало)) составляет всегда * на 1010 * больше, чем pi / 2 перед вашим "средним" углом (вы начинаете с последнего угла перед интервалом) в то время как ваше последнее измерение (Значение (отделка + i-1)) всегда больше , чем pi / 2 позади середины. Таким образом, ваш эффективный диапазон, по которому вы усредняете значение, всегда больше, чем pi, даже если доступны значения данных точно в середине pi / 2 и mid + pi / 2 ... это действительно так?

Так что, если вы действительно заинтересованы в усреднении только значений, где угол меньше, чем pi / 2 от середины, вот мой совет по коду, который, к сожалению, имеет лишь незначительно меньшее время выполнения, чем то, что вы используете в настоящее время. Обратите внимание, что это НЕ рефакторинг, потому что он действует не так, как ваш код, двумя способами, описанными выше.

UnwrappedRadians = unwrap( Angle ./ 180 * pi);
AverageValue = Value;
avgstart=find( UnwrappedRadians > (UnwrappedRadians(1) + pi/2), 2, 'First');
avgend=find( UnwrappedRadians < (UnwrappedRadians(end) - pi/2), 1, 'Last');
for i=avgstart:avgend
    AverageValue(i) = mean(Value(abs(UnwrappedRadians-UnwrappedRadians(i)) <= pi/2));  % nanmean
end

Answer 2

Незначительный рефакторинг сохранит вторую находку в тех случаях, когда вы не найдете результатов, а предварительное распределение с AverageValue со значением сохраняет остальную часть.

UnwrappedRadians = unwrap( Angle ./ 180 * pi);
AverageValue = Value;

for i=1:length(UnwrappedRadians)
    mid = UnwrappedRadians(i);
    start = find( UnwrappedRadians(1:i) < (mid - pi/2), 1, 'Last');
    if ~isempty(start)
        finish = find( UnwrappedRadians(i:end) > (mid + pi/2), 1, 'First');
        if ~isempty(finish)
            AverageValue(i) = mean(Value(start:finish+i-1));  % nanmean
        end
    end
end

Если вы обнаружите, что вычисление финиша пусто чаще, чем вычисление старта, вы можете переключить их порядок так, чтобы проверка финиша выполнялась первой.

Не ясно, будут ли UnwrappedRadians всегда сортироваться. Если это так, вы можете повторно использовать результаты предыдущих поисков, чтобы уменьшить размер подвектора, который вы ищете. Например, если вы ищете последнюю 1 между 11 и концом вектора, а это элемент 23, при поиске 1.1 вы можете уменьшить поиск до 24 и до конца вектора. Я обнаружил, что эта техника может привести к очень большому увеличению скорости.

Векторизация трудна в подобных случаях, потому что вы снова используете индексную переменную (i) в качестве индекса (в операторах поиска). Можно скомпилировать что-то с помощью arrayfun, но вряд ли это будет быстрее (я бы предположил, что медленнее) и определенно будет менее читабельным, чем у вас.