How to resample time vector to regular sampling rate maintaining original values

Hi all,

I am woking with pupil size data collected at 60Hz.

When looking at the derivative of the time vector, I realized that some value have been dropped.

I need to rebuild a timevector with the original time points regularly resampled to 60Hz.

Small violation from 60Hz are fine and I would prefer to keep the original raw timeseries as much as possible.

I am trying to understand how to do this. I tried with some simulated data with griddedInterpolant() but it does not work well.

ANother solution I tried is to take the derivative and whenever it was higher than a threshold, I would fill it with new values consistent with sampling rate. The script works but I cannot believe there is no better way to do this.

Here the code for the second method

time_raw_in_s = Pupil.system_time_stamp/number_microseconds_in_a_second;
sampling_frequency = 60
threshold = round(1/sampling_frequency,3);  
time_full = time_raw_in_s(1);  
for data_point = 1: length(time_raw_in_s)-1
    % disp(['this is trial:', string(data_point)])
    gap = [];
    gap = time_raw_in_s(data_point+1) - time_raw_in_s(data_point); % diff between consecutive time points
    
    if gap > threshold %gap higher than treshold
        % Generate equally spaced points from time_raw_in_s(data_point) to time_raw_in_s(data_point+1)
        new_points = [];
        new_points = time_raw_in_s(data_point):threshold:time_raw_in_s(data_point+1);
        % --- Endpoint check  ---
        % Calculate difference between the last generated point and the actual end
        last_diff =  [];
        last_diff = time_raw_in_s(data_point+1) - new_points(end);
        
        % Handle rounding and near-end edge cases:
        % Sometimes, due to floating-point precision or slightly irregular spacing,
        % the last point may be slightly before or after the intended endpoint.
        if last_diff > 1e-12 && last_diff < threshold/2 % Close enough to end:
            new_points(end) = time_raw_in_s(data_point+1);    %  replace last point with exact endpoint
        elseif last_diff >= threshold/2 % Big enough gap
            new_points = [new_points, time_raw_in_s(data_point+1)]; % append endpoint manually
        end
        % --- Update the complete time vector ---
        % Skip the first point (already included in 'time_full'),
        % and append the new interpolated time points.
        data = [];
        data = new_points(2:end);
        time_full = [time_full; data'];
    elseif gap <= threshold
        % gap < threshold: just add the next original time point
        time_full = [time_full; time_raw_in_s(data_point+1)]
    end 

Answer 1

[David Goodmanson]

Hi Martina,

As I understand it there is a fairly constant clock rate with some values dropped (as opposed to a gap, of say, 2.5 cycles in time that is not commensurate with the clock rate of 60 Hz). If it's the values-dropped case, maybe something like the following:

The vector tnew is the interpolated result, which is displayed below vs the vector of raw data with zeros in the data gaps. All raw data values are preserved. The code probably doesn't work if there are a couple of missing raw data points at either end.

% create raw data with missing points + some random jitter
delt = 1/60;
traw = delt*[1 2 3 5 6 7 9 13 14 15 16 19 20];
traw = traw + .1*delt*rand(size(traw));

% interpolated array
%
% create a longer array of zeros, fill in raw data values
nraw = round(traw/delt);   % indices of raw data
tgaps = zeros(1,nraw(end));
tgaps(nraw) = traw;        % tgaps has value zero for missing data

% find initial and final indices of gaps
tnew = tgaps;
d = [1 diff(nraw)];
fin = nraw(d>1); % index of upper end
wid = d(d>1);    % gap width
sta = fin - wid; % index of lower end
% interpolate
for k = 1:length(fin)
    f = fin(k)
    s = sta(k)
    w = wid(k);
    tnew(s:f) = linspace(tnew(s),tnew(f),w+1);
end

% tnew is the result
format compact
comparison = [tnew;tgaps]
format

comparison =
  Columns 1 through 6
    0.0178    0.0336    0.0502    0.0672    0.0842    0.1016
    0.0178    0.0336    0.0502         0    0.0842    0.1016
  Columns 7 through 12
    0.1172    0.1341    0.1510    0.1675    0.1840    0.2005
    0.1172         0    0.1510         0         0         0
  Columns 13 through 18
    0.2170    0.2346    0.2504    0.2675    0.2843    0.3011
    0.2170    0.2346    0.2504    0.2675         0         0
  Columns 19 through 20
    0.3178    0.3348
    0.3178    0.3348

[Martina]

Hi David,

Thank you so much for this.

I see your point about constant clock.

is there something i can check to see if the asumption applies to my data? Here a picture of two cases of sampling rate drop.

martina

[Martina]

Hi Stephen, Thanks for your feedback. The plot represents the derivative of the time vector. To make it visually clearer I calculated 1/derivative which gives the sampling rate at each time point. Overall it is around 60 hz but whenever it goes down it means the gap between two consecutive points is larger than 0.016s which corresponds to a sampling rate of 60 more or less. Hope the graph is clearer now Thanks so much for ur comments Martina