大家好，又见面了，我是你们的朋友全栈君。
0, 说明
关于插值，官网有个小总结，可以直接去参考（从1维到多维），下面是我举的例子。
1, 一维插值interp1(x,y,X1,method)
x = linspace(0,10,11)
y = sin(x)
plot(x,y,'-ro')
插值方法有如下：
method=‘nearest’,‘linear’,‘spline’,‘pchip’,‘cubic’
比如使用三次条样插值spline，则
x = linspace(0,10,11)
y = sin(x)
plot(x,y,'-ro')
xnew = linspace(0,10,101)
f = interp1(x,y,xnew,'spline')
plot(xnew,f)
2,高维插值
2.1 二维插值
使用interp2(X,Y,Z,X1,Y1,method)函数, method包含：nearest,‘linear’,‘spline’,‘cubic’
x = linspace(-5,5,11);
y = linspace(-5,5,11);
[X,Y] = meshgrid(x,y);
Z = 3.*sin(X.^2+Y.^2);
xi = linspace(-5,5,100);
yi = linspace(-5,5,100);
[XI,YI]=meshgrid(xi,yi);
ZI = interp2(X,Y,Z,XI,YI,'cubic');
surf(XI,YI,ZI),view(-25,25)
2.2 三维插值
参考官网，用法如下
Vq = interp3(X,Y,Z,V,Xq,Yq,Zq)
Vq = interp3(V,Xq,Yq,Zq)
但这个方法有时候不是很实用，主要是X,Y,Z都要弄成网格的形式，因此我一般大于等于三维的时候就直接采用interpn了，见下一节。
2.3 多维插值
大于等于三维一般采用interpn,可以参考官网用法，
Vq = interpn(X1,X2,...,Xn,V,Xq1,Xq2,...,Xqn)
这里x1,x2,...xn,以及xq1,xq2,...xqn都是一维的（注意是单调函数），所以比Interp3简单很多。
方法有如下method: ‘linear’, ‘nearest’, ‘pchip’,‘cubic’, ‘makima’, or ‘spline’. The default method is ‘linear’.
举例：
1)插值一个点
现在有一个高维数据（4维），横坐标是经度，纵坐标是维度，高是海拔，V的值是在这三维中的水汽含量。我现在有了V的数据，这个数据是(37,10,10)的大小，表示高有37层，经纬度分别都是10的大小（因此经纬度构成100的数据网格），现在要计算高500m，经纬度分别为（80，32）的点的值（插值）
data_path = 'C:\Users\Lenovo\Desktop\recent work\matlab_interp\'
el = importdata([data_path,'el_sort.txt']); %海拔
lat = importdata([data_path,'lat_sort.txt']); %纬度
lon = importdata([data_path,'lon_sort.txt']); %经度
V = readNPY([data_path,'s_h_reverse.npy']); %数据
Vq = interpn(el,lon,lat,V,500,80,32,'spline') %插值
2）插值两个点
上面插值只在一个点（500，80，30）上进行，但有时我们要插值的是很多个点构成的数组。
%在(500,80.03,32.3)和(900,81.2,31.4)两个点插值
Vq = interpn(el,lon,lat,V,[500 900],[80.03 81.2],[32.3 31.4],'spline')
- 多个点（数组）
Vq = interpn(el,lon,lat,squeeze(s_h_full(i,:,:,:)),coor(:,1),coor(:,2),coor(:,3),'spline')
这里coor表示坐标，有海拔，经，纬度，都是一列。
完整代码
data_path = 'C:\Users\Lenovo\Desktop\recent work\matlab_interp\';
el = importdata([data_path,'el_sort.txt']); %海拔
lat = importdata([data_path,'lat_sort.txt']); %纬度
lon = importdata([data_path,'lon_sort.txt']); %经度
V = readNPY([data_path,'s_h_full_inv.npy']); %数据
coor = readNPY([data_path,'coor_in_new.npy']);
levell = 50;
hours = 720;
s_h_interp = [];
for i = 1:hours
Vq = interpn(el,lon,lat,squeeze(V(i,:,:,:)),coor(:,1),coor(:,2),coor(:,3),'spline');
s_h_interp = horzcat(s_h_interp,Vq);
end
writeNPY(s_h_interp,[data_path,'s_h_interp_matlab_161.npy'])
发布者：全栈程序员-站长，转载请注明出处：https://javaforall.net/136155.html原文链接：https://javaforall.net
