# Plotting

figure(1);

scatter(rocket.t,rocket.v,100,’fill’);

grid(‘minor’):

xlabel(‘t’);

ylabel(‘v’);

# Contents

Detailed Interpolation of a Single Point

The Interpolation Functions interp1

# Detailed Interpolation of a Single Point

Suppose we have a rocket and have recorded it’s speed during the list of times.

And we want to know it’s speed at time t=16 s. There are a number of ways we can estimate this value or interpolate this value from the given data.

## Nearest Point Interpolation (1 data point)

Simply by looking at the data we can see that the nearest time point to t=16 s is t=15 s.

We can take the speed at t=15 s and use it as a first hand estimate for the speed at t=16 s and this is known as single point interpolation.

`St16=St15`

St16=362.78

## Linear Interpolation (2 data points)

Suppose we want to be a bit more accurate we can instead take two points and calculate the equation for a line. We can then plug in our x data point and use the equation of the straight line to estimate our unknown y data point. Looking at the data for t=16 s the nearest two time points are t=15 s and t=20 s.

The equation for a straight line is:

For two times we can write this out in matrix form:

Now we can substitute in the known values:

We can save the known data as variables S and T:

`S=[362.78;517.35]`

`T=[1,15;1,20]`

And call the unknown variable U:

We know that:

However we want to solve for U. We can get U on it’s own by multiplying throughout by the inverse matrix of T:

Since a matrix multiplied by it’s inverse gives the identity matrix and a matrix multipied by the identity matrix gets unchanged:

Which can be reordered:

Thus to solve for U we need the inverse matrix of T:

This can be calculated using the function `Minv=inv(M)`

`Tinv=inv(T)`

Thus we can get the unknown coefficients of U using:

`U=Tinv*S`

i.e. in the limits of 15<t<20 a0=-100.93 and a1=30.914

Instead of using the inverse matrix we could also use left division

`U=T\S`

Now that we have the coefficients we can calculate the S(16)

`St16=[1,16]*U`

We can check the equation using the known data:

Unsurprisingly the values at the limits t=15 and t=20 give exact matches. There is a slight deviation just outside this interval and a massive deviation substantially away from this interval.

## Quadratic Interpolation (3 data points)

Suppose we want to be a bit more accurate we can instead take three points and calculate a quadratic equation. We can then plug in our x data point and use the quadratic equation to estimate our unknown y data point. Looking at the data for t=16 s the nearest three time points are t= 10 s, t=15 s and t=20 s.

The equation for a quadratic function is:

For three times we can write this out in matrix form:

Now we can substitute in the known values:

We can save the known data as variables S and T:

`S=[227.04;362.78;517.35]`

`T=[1,10,10^2;1,15,15^2;1,20,20^2]`

And call the unknown variable U:

To solve for U we use:

`U=T\S`

i.e. in the limits of 10<t<20 a0=12.05, a1=17.733 and a2=0.3766

Now that we have the coefficients we can calculate the S(16)

`St16=[1,16,16^2]*U`

## Cubic Interpolation (4 data points)

Suppose we want to be a bit more accurate we can instead take four points and calculate a cubic equation. We can then plug in our x data point and use the cubic equation to estimate our unknown y data point. Looking at the data for t=16 s the nearest three time points are t= 10 s, t=15 s, t=20 s and t=22.5 s.

The equation for a cubic function is:

For four times we can write this out in matrix form:

Now we can substitute in the known values:

We can save the known data as variables S and T:

`S=[227.04;362.78;517.35;602.97]`

`T=[1,10,10^2,10^3;1,15,15^2,15^3;1,20,20^2,20^3;1,22.5,22.5^2,22.5^3]`

And call the unknown variable U:

To solve for U we use:

`U=T\S`

i.e. in the limits of 10<t<22.5 a0=4.254, a1=21.2655333, a2=0.13204 and a3=0.0054347

Now that we have the coefficients we can calculate the S(16)

`St16=[1,16,16^2,16^3]*U`

## Comparison

We can compare our interpolation for t=16 using the four different methods. As we can see when we are interpolating with higher order polynomials the result begins to converge as expected:

# Interpolation of a Data Set

As we seen it is quite involved interpolating the data for a single data point. If we want to interpolate 100’s or even 1000’s of data points manually going through them all like above it would take all day. Luckily Octave/MATLAB have inbuilt interpolation functions.

`intery=interp1(originalx,originaly,xdatapoints,'Method')`

We can type in the data using Octave/MATLAB.

[code language=”matlab”]

data=[0,0;

10,227.04;

15,362.78;

20,517.35;

22.5,602.97;

30,901.67]

[/code]

We could also use Excel to copy and paste in the data. In MATLAB which is more polished, we can directly go to the variable editor, right click it and create a new variable and then name accordingly. This option isn’t available in Octave however we can quickly create a variable, assign in to a scalar and then copy and paste in a matrix in the variable editor.

For instance typing in:

`data=1;`

Then clicking the new variable data in the workspace to open it up in the variable editor:

We can right click the top cell of data in the variable editor and we can right click Paste Table to copy the Excel Data through to Octave:

Now we have our data variable:

Next we need to make the new x data. For this I am going to make a new variable called dataint and I want the first column to go from 1 to 30. Here logical indexing is used when creating the variable dataint. The `(:,1)`

means we are creating data in all rows but using the first column. The `[1:30]'`

is a column vector ranging from 1 to 30 in steps of 1. If you are unsure about logical indexing see my earlier guide Introduction to Arrays.

`dataint(:,1)=[1:30]'`

Now we are looking to interpolate the data:

`intery=interp1(originalx,originaly,xdatapoints,'Method')`

Our original x data originalx is the first column of data our original y data originaly is the second column of data

`originalx=data(:,1)`

`originaly=data(:,2)`

Our x datapoints `xdatapoints`

are the first column in dataint

`xdatapoints=dataint(:,1)`

We will start with Nearest Interpolation meaning instead of typing in `'Method'`

Â we type in `'Nearest'`

. We can now get intery:

`intery=interp1(originalx,originaly,xdatapoints,'Nearest')`

We can clear up our code here:

`intery=interp1(originalx,originaly,xdatapoints,'Nearest')`

By substituting these in directly:

`originalx=data(:,2)`

`originaly=data(:,2)`

`xdatapoints=dataint(:,1)`

This becomes

`intery=interp1(data(:,1),data(:,2),dataint(:,1),'Nearest')`

Now it is more likely we want the y data to be alongside the xdatapoints so instead of having the separate output variable `intery`

we can add it as a second column to `dataint`

Â by using:

`dataint(:,2)=interp1(data(:,1),data(:,2),dataint(:,1),'Nearest')`

We can press the [â†‘] key to get our last line of code and easily modify the line of code to print linear interpolated data to column 3. For linear interpolation the method is `'Linear'`

:

`dataint(:,3)=interp1(data(:,1),data(:,2),dataint(:,1),'Linear')`

MATLAB/Octave’s `interp1`

function doesn’t have an inbuilt quadratic interpolation so I will go directly to cubic interpolation and here the method is `'Pchip'`

Â , this time I want the data in the 4th column.

`dataint(:,4)=interp1(data(:,1),data(:,2),dataint(:,1),'Pchip')`

Since we were interested at the value at t=16 which is the 16th row we can output the 16th row into the command window:

`t16row=`

`dataint(16,:)`

These are very similar to the values manually calculated, the differences are likely due to rounding.

The code above can be cleared up and put in a script file.

[code language=”matlab”]

% Data

data=[0,0;

10,227.04;

15,362.78;

20,517.35;

22.5,602.97;

30,901.67];

% Lower and Upper Limits for interpolation

m=1;

n=30;

% Interpolation Spreadsheet

dataint(:,1)=[m:n]’;

dataint(:,2)=interp1(data(:,1),data(:,2),dataint(:,1),’Nearest’);

dataint(:,3)=interp1(data(:,1),data(:,2),dataint(:,1),’Linear’);

dataint(:,4)=interp1(data(:,1),data(:,2),dataint(:,1),’Pchip’);

[/code]