for loop jump an element of an array
23 Ansichten (letzte 30 Tage)
Ältere Kommentare anzeigen
Hello there,
it's my first post here, I'm writing a code to perform linear interpolation of a sinusoid and then a resampling with an higeher sampling frequency but i'm facing some trouble since the foor loop seems to jump one value of the sampling time vector and that's a problem when i'm walking through the array to find the overlapping times of the two sampling time vectors.
f0 = 10; % Hz
T0 = 1/f0; % period of sin [s]
fc = 100; %Hz
T = 1/fc; % first sampling period [s]
N = floor(10*(fc/f0)); % number of samples in 10 periods
samples_time = 0:T:(N-1)*T;
samples_value = zeros(1,N);
% first sampling
for i=1:N
samples_value(i) = sin(2*pi*f0*samples_time(i));
end
fc2 = 200; % Hz
T2 = 1/fc2; % second sampling period [s]
N2 = floor((N-1)*(T/T2)+1); % number of samples
samples_time2 = 0:T2:(N2-1)*T2;
for j=1:N2
if samples_time2(j) == 0.64
disp("found");
end
%disp(samples_time2(j));
end
The foor loop seems not to find the value 0.64 in the array, it works with other values. I cannot find a solution, it's very strange. (I'm running R2023b)
Hope someone can help me, thanks
1 Kommentar
Stephen23
am 27 Mär. 2024
"I cannot find a solution, it's very strange."
It is not strange at all: you made the mistake of testing for exact equality of binary floating point numbers.
The solution is also very simple: always compare the absolute difference against a tolerance:
abs(A-B)<tol
More information on this topic:
This is worth reading as well:
Antworten (2)
Voss
am 27 Mär. 2024
"x = j:i:k creates a regularly-spaced vector x using i as the increment between elements. The vector elements are roughly equal to [j,j+i,j+2*i,...,j+m*i] where m = fix((k-j)/i). However, if i is not an integer, then floating point arithmetic plays a role in determining whether colon includes the endpoint k in the vector, since k might not be exactly equal to j+m*i."
I'd expand that to add, if i is not an integer, then floating point arithmetic plays a role in determining whether colon includes any intended value beyond the first.
To avoid this problem, construct your time vectors using a spacing of 1 in colon, i.e.:
samples_time = (0:N-1)*T;
samples_time2 = (0:N2-1)*T2
instead of:
samples_time = 0:T:(N-1)*T;
samples_time2 = 0:T2:(N2-1)*T2;
(Yes, the same problem of missing values happens in samples_time as well.)
Some illustrations:
f0 = 10; % Hz
T0 = 1/f0; % period of sin [s]
fc = 100; %Hz
T = 1/fc; % first sampling period [s]
N = floor(10*(fc/f0)); % number of samples in 10 periods
Constructing samples_time both ways:
samples_time = 0:T:(N-1)*T
samples_time_test = (0:N-1)*T
They look the same, but they are not:
isequal(samples_time,samples_time_test)
In particular, 0.64 is in the one constructed using a colon spacing of 1 but not the one using T spacing:
find(samples_time == 0.64)
find(samples_time_test == 0.64)
There is an element in samples_time that is very close to 0.64 (about 1e-16 away), but it's not exactly 0.64:
[d,idx] = min(abs(samples_time-0.64))
whereas in samples_time_test that element is exactly 0.64:
[d,idx] = min(abs(samples_time_test-0.64))
The same is true for the second set of times.
fc2 = 200; % Hz
T2 = 1/fc2; % second sampling period [s]
N2 = floor((N-1)*(T/T2)+1); % number of samples
samples_time2 = 0:T2:(N2-1)*T2
samples_time2_test = (0:N2-1)*T2
isequal(samples_time2,samples_time2_test)
find(samples_time2 == 0.64)
find(samples_time2_test == 0.64)
And 0.64 isn't the only missing value either, in either vector:
% elements in samples_time_test that aren't in samples_time:
samples_time_test(~ismember(samples_time_test,samples_time))
% elements in samples_time2_test that aren't in samples_time2:
samples_time2_test(~ismember(samples_time2_test,samples_time2))
4 Kommentare
Stephen23
am 29 Mär. 2024 um 7:59
"...can I have the same problem?"
Yes. Using two different methods can result in two different values. That is nature of binary floating point numbers.
the cyclist
am 27 Mär. 2024
Instead, check equality with a small tolerance, e.g.
if abs(samples_time2(j) - 0.64) < 1.e-6
0 Kommentare
Siehe auch
Kategorien
Mehr zu Logical finden Sie in Help Center und File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!