Efficiently calculating sum-thresholds across vector

3 views (last 30 days)
Joey Costello on 22 Sep 2020
Edited: Bruno Luong on 22 Sep 2020
I have a 1xN vector of samples, for which I need to calculate the indices of when a moving sum exceeds a threshold.
One could think of the samples as current that is charging a capacitor, and once the capacitor is charged to the threshold, a pulse is fired and the capacitor is reset. I need to find the times when the pulse is fired.
I have a working loop implementation below, however, it is relatively slow and I'd like to find a more efficient method. Is there a vectorized or built-in method for finding these pulse times?
samples;
threshold = 50;
pulseIdx = zeros(size(samples));
runningSum = 0;
for i = 1:length(samples)
runningSum = runningSum + samples(i);
if runningSum > threshold
pulseIdx(i) = 1;
runningSum = runningSum - threshold; % allow spillover to next cycle
end
end
Turlough Hughes on 22 Sep 2020
Are the values in samples always positive?

Turlough Hughes on 22 Sep 2020
Edited: Turlough Hughes on 22 Sep 2020
You could do the following but with the caveat that it assumes samples are all positive values
runningSum = cumsum(samples);
idx = find(diff(mod(runningSum,threshold))<0)+1;
Or in the case where they are not all positive, perhaps something like the following:
idx = find(diff(mod(runningSum,50))<-0.5*threshold)+1

Ameer Hamza on 22 Sep 2020
Edited: Ameer Hamza on 22 Sep 2020
Try something like this
rng(0);
samples = rand(1, 1000);
threshold = 10;
samples_sum = cumsum(samples);
idx = [];
while any(samples_sum > threshold)
i = find(samples_sum > threshold, 1);
idx = [idx i];
samples_sum = samples_sum - samples_sum(i-1);
end
idx are indexes of samples where sum exceed threshold.

Bruno Luong on 22 Sep 2020
Edited: Bruno Luong on 22 Sep 2020
This works regardless the sign of sample
% Test data
samples = 5*rand(1,10000);
samples = samples+0.5*randn(size(samples)); % with noise
threshold = 50;
pulseIdx = zeros(size(samples));
runningSum = 0;
for i = 1:length(samples)
runningSum = runningSum + samples(i);
if runningSum > threshold
pulseIdx(i) = 1;
runningSum = runningSum - threshold; % allow spillover to next cycle
end
end
% My method
cs = cummax(cumsum(samples));
n = floor(cs(end)/threshold);
bdr = (1:n)*threshold;
[~,~,loc]=histcounts(bdr,[0,cs]);
pulseIdx2 = zeros(size(samples));
pulseIdx2(loc(loc>0)) = 1;
% This should return 1, excepted when numerical errors spoil thing
isequal(pulseIdx,pulseIdx2)

R2018a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!