for loop jump an element of an array

28 views (last 30 days)
Davide
Davide on 27 Mar 2024
Commented: Davide on 30 Mar 2024 at 11:29
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

Sign in to answer this question.

Answers (2)

Voss
Voss on 27 Mar 2024
Ran in:
From the documentation for the colon operator:
"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 = 1×100
0 0.0100 0.0200 0.0300 0.0400 0.0500 0.0600 0.0700 0.0800 0.0900 0.1000 0.1100 0.1200 0.1300 0.1400 0.1500 0.1600 0.1700 0.1800 0.1900 0.2000 0.2100 0.2200 0.2300 0.2400 0.2500 0.2600 0.2700 0.2800 0.2900
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
samples_time_test = (0:N-1)*T
samples_time_test = 1×100
0 0.0100 0.0200 0.0300 0.0400 0.0500 0.0600 0.0700 0.0800 0.0900 0.1000 0.1100 0.1200 0.1300 0.1400 0.1500 0.1600 0.1700 0.1800 0.1900 0.2000 0.2100 0.2200 0.2300 0.2400 0.2500 0.2600 0.2700 0.2800 0.2900
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
They look the same, but they are not:
isequal(samples_time,samples_time_test)
ans = logical
0
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)
ans = 1×0 empty double row vector
find(samples_time_test == 0.64)
ans = 65
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))
d = 1.1102e-16
idx = 65
whereas in samples_time_test that element is exactly 0.64:
[d,idx] = min(abs(samples_time_test-0.64))
d = 0
idx = 65
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 = 1×199
0 0.0050 0.0100 0.0150 0.0200 0.0250 0.0300 0.0350 0.0400 0.0450 0.0500 0.0550 0.0600 0.0650 0.0700 0.0750 0.0800 0.0850 0.0900 0.0950 0.1000 0.1050 0.1100 0.1150 0.1200 0.1250 0.1300 0.1350 0.1400 0.1450
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
samples_time2_test = (0:N2-1)*T2
samples_time2_test = 1×199
0 0.0050 0.0100 0.0150 0.0200 0.0250 0.0300 0.0350 0.0400 0.0450 0.0500 0.0550 0.0600 0.0650 0.0700 0.0750 0.0800 0.0850 0.0900 0.0950 0.1000 0.1050 0.1100 0.1150 0.1200 0.1250 0.1300 0.1350 0.1400 0.1450
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
isequal(samples_time2,samples_time2_test)
ans = logical
0
find(samples_time2 == 0.64)
ans = 1×0 empty double row vector
find(samples_time2_test == 0.64)
ans = 129
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))
ans = 1×13
0.6400 0.6500 0.6600 0.6700 0.6800 0.6900 0.7000 0.8200 0.8300 0.9200 0.9300 0.9400 0.9500
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% elements in samples_time2_test that aren't in samples_time2:
samples_time2_test(~ismember(samples_time2_test,samples_time2))
ans = 1×28
0.5150 0.5250 0.5350 0.5450 0.5550 0.5650 0.5750 0.6400 0.6500 0.6600 0.6700 0.6800 0.6900 0.7000 0.7850 0.7950 0.8050 0.8150 0.8200 0.8250 0.8300 0.9200 0.9300 0.9350 0.9400 0.9450 0.9500 0.9550
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
  4 Comments
Show 2 older comments
Stephen23
Stephen23 on 29 Mar 2024 at 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.
Davide
Davide on 30 Mar 2024 at 11:29
Sure, of course. As soon as I get home I will provide you the complete code

Sign in to comment.

the cyclist
the cyclist on 27 Mar 2024
You should not use strict equality to compare floating point numbers.
Instead, check equality with a small tolerance, e.g.
if abs(samples_time2(j) - 0.64) < 1.e-6

Categories

Find more on Logical in Help Center and File Exchange

Tags

Community Treasure Hunt

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

Start Hunting!