Efficient indexing of bounded matrices without for loop

Hi there,

I have a matrix of tide data at irregular sub-minute intervals that includes time and water level heights (gMedt). Time is in the first column and water level in the second. I'm basically trying to resample the data to hourly data (tA1 and tA2) from the first available hour to the last available hour in two different ways as follows. Columns 3 and 4 represents the finer resolution time and water level just before the hourly mark and columns 5 and 6 represent the same just after the hourly mark. Column 2 would then be a linear interpolation between the two water levels of Columns 4 and 6.

res = 1/24;
tAst = ceil( gMedt( 1, 1 ) / res ) * res;
tAen = floor( gMedt( end, 1 ) / res ) * res;
tA1 = tAst:res:tAen; tA1 = tA1';
tA2 = tA1;
for i = 1:length( tA1 );
  tA1( i, 3 ) = gMedt( numel( gMedt( gMedt( :, 1 ) < tA1( i, 1 ) ) ), 1 );
  tA1( i, 4 ) = gMedt( numel( gMedt( gMedt( :, 1 ) < tA1( i, 1 ) ) ), 2 );
  tA1( i, 5 ) = gMedt( 1+numel( gMedt( gMedt( :, 1 ) < tA1( i, 1 ) ) ), 1 );
  tA1( i, 6 ) = gMedt( 1+numel( gMedt( gMedt( :, 1 ) < tA1( i, 1 ) ) ), 2 );

tA2( i, 3 ) = gMedt( find( gMedt( :, 1 ) < tA2( i, 1 ), ...
    1, 'last' ), 1 );
tA2( i, 4 ) = gMedt( find( gMedt( :, 1 ) < tA2( i, 1 ), ...
    1, 'last' ), 2 );
tA2( i, 5 ) = gMedt( find( gMedt( :, 1 ) > tA2( i, 1 ), ...
    1, 'first' ), 1 );
tA2( i, 6 ) = gMedt( find( gMedt( :, 1 ) > tA2( i, 1 ), ...
    1, 'first' ), 2 );
end

Neither of these ways is particularly fast especially since the gMedt variable is 16 million rows long and tA1 and tA2 are 43825 rows long. Is there a way to index the tA matrix by finding the correct gMedt row without a for loop?

Answer 1

[Walter Roberson]

Do not repeat the calculations, and do not index when you do not need to.

For example,

mask = 
posn = sum(mask);    %same as numel( gMedt(mask) )

tA1( i, 3 ) = gMedt( posn, 1 );
tA1( i, 4 ) = gMedt( posn, 2 );
tA1( i, 5 ) = gMedt( 1+posn, 1 );
tA1( i, 6 ) = gMedt( 1+posn, 2 );

This can in turn be shortened to

posn = sum( gMedt( :, 1 ) < tA1( i, 1 ) );
tA1(i, [3 4]) = gMedt( posn, 1:2 );
tA1(i, [5 6]) = gMedt( 1 + posn, 1:2 );

[curoi]

Walter,

It looks like the resample function might actually be quicker. However, it excludes NaNs from the interpolation which is not what I want. Is there resample type function that doesn't exclude NaN values. For example, it will only interpolate to resampled times if that time doesn't fall between an original recording and NaN value. If the resampled time does fall between an original time and a NaN value, then it becomes NaN.

[Walter Roberson]

interp1() and interp2() will cheerfully produce nans.

[curoi]

I don't just want it to produce NaN's, I'd like it to accept NaNs in the input. I've checked the file exchange and there are some functions that will interpolate over NaN inputs. However, what I'd like is the input NaNs to be used as NaN producers. If the output time falls between a NaN input and a non-NaN input time, the water level output should be NaN.