Speeding up interpolation with scattered interpolation

Hi,

I've been working on speeding up my code for a project and I've identified that I'm spending a lot of time interpolating my data. I have data for wind speeds, air temperature, pressure and a few other things that depend on position in the X,Y and Z directions. I've been using ScatteredInterpolation to interpolate across 1425 data points and then using the interpolant find the wind speeds, air temperature, pressure etc. at a specific position. I was wondering if there was a faster way to do this, either speeding up scatteredinterpolation or a different interpolating method. Am fairly new to MATLAB so help appreciated.

Thanks.

[daniel]

I've just noticed that actually using the interpolant to find my wind speed, air tempearture etc. at the desired location is also taking up a lot of time (60% of my overall run time) whilst the scattered interpolation takes up around 35% of my overall run time.

[Matt J]

We need to be able to run your code to say anything.

[daniel]

The whole code is very long, slighly messy and would be a pain to share on here. If its possible not to share the whole thing I'd rather do that if thats okay. Ive attached below the sections of code for the interpolation.

The Atmospheric data arrays contain data from the GFS weather model, more specifically they contain Temperature, Pressure, Air Density, Eastwards windspeed, northwards windspeed and vertical windspeed at a range of points.

The points are spaced every 0.25 lat and long and there is data for 57 different altitudes.

When I am interpolating I'm interpolating across a 5 latitudes and 5 longitudes and all altitudes.

AtmosphericDataPrevious = G.AtmosphericDataPrevious;
AtmosphericDataCurrent = G.AtmosphericDataCurrent;
AtmosphericDataNext = G.AtmosphericDataNext;
PreviousForecastTime = str2double(G.PreviousForecastTime);
CurrentForecastTime = str2double(G.CurrentForecastTime);
NextForecastTime = str2double(G.NextForecastTime);
AltitudeData = reshape(AtmosphericDataCurrent(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,1),[1425,1]);
LatData0 = LowerLat:0.25:UpperLat;
LatData1 = repmat(LatData0,[5,1,57]);
LatData = reshape(LatData1,[1425,1]);
LonData0 = LeftLon:RightLon;
LonData1 = reshape(LonData0,[5,1]);
LonData2 = repmat(LonData1,[1,5,57]);
LonData = reshape(LonData2,[1425,1]);
%previous data
TemperatureArrayPrevious = reshape(AtmosphericDataPrevious(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,4),[1425,1]);
TemperaturePrevious0 = scatteredInterpolant(LatData,LonData,AltitudeData,TemperatureArrayPrevious,"natural");
DensityArrayPrevious = reshape(AtmosphericDataPrevious(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,3),[1425,1]);
DensityPrevious0 = scatteredInterpolant(LatData,LonData,AltitudeData,DensityArrayPrevious,"natural");
PressureArrayPrevious = reshape(AtmosphericDataPrevious(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,2),[1425,1]);
PressurePrevious0 = scatteredInterpolant(LatData,LonData,AltitudeData,PressureArrayPrevious,"natural");
VerticalWindVelocityArrayPrevious = reshape(AtmosphericDataPrevious(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,5),[1425,1]);
VerticalWindVelocityPrevious0 = scatteredInterpolant(LatData,LonData,AltitudeData,VerticalWindVelocityArrayPrevious,"natural");
EastwardsWindVelocityArrayPrevious = reshape(AtmosphericDataPrevious(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,6),[1425,1]);
EastwardsWindVelocityPrevious0 = scatteredInterpolant(LatData,LonData,AltitudeData,EastwardsWindVelocityArrayPrevious,"natural");
NorthwardsWindVelocityArrayPrevious = reshape(AtmosphericDataPrevious(LowerLatLocation:UpperLatLocation,LeftLonLocation:RightLonLocation,:,7),[1425,1]);
NorthwardsWindVelocityPrevious0 = scatteredInterpolant(LatData,LonData,AltitudeData,NorthwardsWindVelocityArrayPrevious,"natural");
TemperaturePrevious = TemperaturePrevious0(CurrentLatitude,CurrentLongitude,CurrentAltitude);%Temperature at altitude (K)
DensityPrevious = DensityPrevious0(CurrentLatitude,CurrentLongitude,CurrentAltitude);%Air density at altitude (Kg/m^3)
PressurePrevious = PressurePrevious0(CurrentLatitude,CurrentLongitude,CurrentAltitude); %Pressure at altitude (Pa) 
VerticalWindVelocityPrevious = VerticalWindVelocityPrevious0(CurrentLatitude,CurrentLongitude,CurrentAltitude); %Vertical wind-velocity (m/s)
EastwardsWindVelocityPrevious = EastwardsWindVelocityPrevious0(CurrentLatitude,CurrentLongitude,CurrentAltitude);%Eastwards wind-velocity (m/s)
NorthwardsWindVelocityPrevious = NorthwardsWindVelocityPrevious0(CurrentLatitude,CurrentLongitude,CurrentAltitude); %Northwards wind-velocity (m/s)

Once I've interpolated im using the interpolant to find the atmospheric properties at a given point.

[Torsten]

The points are spaced every 0.25 lat and long and there is data for 57 different altitudes.

If this is the case, why don't you use GriddedInterpolant instead of ScatteredInterpolant ? The interpolation afterwards will be much faster.

Answer 1

[Matt J]

Your data is gridded, so you should be using griddedInterpolant, as opposed to scatteredInterpolant. Gridded interpolation is much faster than scattered interpolation.

TemperaturePrevious0 = griddedInterpolant( {LatData0,LonData0,AltitudeData0},   TemperatureArrayPrevious)

Also...

Once I've interpolated im using the interpolant to find the atmospheric properties at a given point.

...you should not interpolate one point at a time. You should supply a vector of query points and interpolate them in a single call.

[Matt J]

As an additional remark, your use of repmat is unnecessarily cumbersome. You could have derived LatData and LonData in a single line using ndgrid,

[LatData,LonData]=ndgrid( LowerLat:0.25:UpperLat   ,  LeftLon:RightLon)

although note that that too would be unnecessary. None of Matlab's interpolation routines actually require the data to be in ndgrid form.

[daniel]

I've tried the above code but i keep getting the following error:

Error using griddedInterpolant

Interpolation requires at least two sample points for each grid dimension.

With this; You should supply a vector of query points and interpolate them in a single call.

I don't think I'm able to supply a vector of query points. The temperature, wind speed etc. are being used to solve an ode using a 4th order runge-kutta method, and my ODE depends on these values. So the next point that I'd want to interpolate is unknow until i have interpolated the previous point and used those values in my ODE. (hopefully i explained that well enough).

[Matt J]

I've tried the above code but i keep getting the following error:

It means one of the coordinates you are interpolating (Lat, Lon, Altitude) has less than two sample points.

[daniel]

I've got it sorted now, thank you very much for your help. sorry if my lack of knowledge was frustrating.