Looping interpolation scheme has errors for some (but not all) datasets?

Hello,

I'm having trouble with an interpolation function I'm trying to write, for use with multiple large datasets. For each dataset, I will have a given time period and two time series: 'a', with data taken at 64-second intervals, and 'b', with data taken at 60-second intervals. I want to resample the 'a' dataset to have the same length as the 'b' dataset, using this code:

resampled_a_data = 1:1:length(b_data);
a_timestamps = 0:64:(64*length(a_data));
timestamp = 0;
for i = resampled_a_data
    %find index j of closest value in a_timestamps 
    j = interp1(a_timestamps,1:length(a_timestamps),timestamp,'nearest');
    a_datapoint = a_data(j);
    %replace resampled_a_data[i] with a_data[j]
    resampled_a_data(i) = a_datapoint;
    timestamp = timestamp + 60;
end

I have tried this code with three datasets so far. One of them works perfectly well, but for the other two datasets, the scheme gives an error at

a_datapoint = a_data(j);

when i is somewhere in the 40000s (it is different for both datasets), and for all i after that datapoint. Prior to giving an error, the contents of the array resampled_a_data are exactly what I would expect them to be. Can anyone advise me on how to fix this error, or alternatively, a better way to accomplish the same goal of resampling the datasets?

Answer 1

[Voss]

Notice that this:

a_timestamps = 0:64:(64*length(a_data))

or equivalently this:

a_timestamps = (0:length(a_data))*64

makes a_timestamps a vector with length length(a_data)+1

For instance:

a_data = 1:100;                         % 1-by-100
numel(a_data)

ans = 100

a_timestamps = (0:length(a_data))*64;   % 1-by-101
numel(a_timestamps)

ans = 101

So changing it to this:

a_timestamps = (0:length(a_data)-1)*64; % 1-by-100
numel(a_timestamps)

ans = 100

might fix the problem (hard to say for sure without having the data or knowing what the error message is).

[Voss]

However, you can likely avoid the loop entirely and do all your interpolations at once:

% some made-up data for demonstration
a_data = cos((1:10)/2);
b_data = cos((1:ceil(numel(a_data)*64/60))/2);

% define the timestamps as above
a_timestamps = (0:length(a_data)-1)*64;
b_timestamps = (0:length(b_data)-1)*60;

% resample (two different ways)
nearest_a_data = interp1(a_timestamps,a_data,b_timestamps,'nearest','extrap');
linear_a_data = interp1(a_timestamps,a_data,b_timestamps,'linear','extrap');

% plot
plot(a_timestamps,a_data,'-o','LineWidth',2)
hold on
plot(b_timestamps,b_data,'-o','LineWidth',2)
plot(b_timestamps,nearest_a_data,'--x','LineWidth',2)
plot(b_timestamps,linear_a_data,'--x','LineWidth',2)
legend({'a','b','a\_resampled (nearest)','a\_resampled (linear)'})

[Chloe Liebenthal]

Replacing my code with this method worked perfectly. Thank you so much for helping me figure this out. I really appreciate it!

[Voss]

You're welcome!