How to get an equation of spline interpolation for a set of data X and Y?

8 次查看（过去 30 天）
zaid jamil
zaid jamil 2020-6-18
评论： Fayyaz Ahmad 2021-10-20
Hi , i have a set of data which is x and y and lets say for example
x=[1 2 3 4 5 6 7 8 9 10 ]
y=[100 200 300 400 500 600 700 800 900 1000]
How to get the equation of spline interpolation for this data, as i used the command
OT = spline(x,y);
  2 个评论
zaid jamil
zaid jamil 2020-6-18
ok , OT is a struct that have many variables wich contains coeffecants and breaks and many other things... how i can use it to make an equation so that by using it in anther program it can work perfectly

请先登录，再进行评论。

请先登录，再回答此问题。

回答（2 个）

John D'Errico
John D'Errico 2020-6-18
编辑：John D'Errico 2020-6-18
Congratulations! You win the door prize as the person to ask this question the millionth time!
Seriously, this is a question I get asked so many times I have lost count. And there really is no easy formula you can write down, because the answer is a bit of a mess. I do offer a code to provide what you want. But reember that a spline is a segmented thing, composed of segments of cubic polynomials, defined in a piecewise sense.
x = linspace(0,2*pi,9)
x =
0 0.7854 1.5708 2.3562 3.1416 3.927 4.7124 5.4978 6.2832
y = sin(x);
spl = spline(x,y);
As I said, spl is a cubic spline. In my SLM toolbox, I provide a code that will extract the polynomial segments.
funs = slmpar(spl,'symabs')
funs =
2×8 cell array
Columns 1 through 7
{1×2 double} {1×2 double} {1×2 double} {1×2 double} {1×2 double} {1×2 double} {1×2 double}
{1×1 sym } {1×1 sym } {1×1 sym } {1×1 sym } {1×1 sym } {1×1 sym } {1×1 sym }
Column 8
{1×2 double}
{1×1 sym }
So in the half open interval [0,pi/4), the cubic polynomial is given as
>> funs{1,1}
ans =
0 0.7854
>> vpa(funs{2,1},10)
ans =
- 0.08493845396*x^3 - 0.1356173501*x^2 + 1.059224243*x
In the half open interval [pi, pi + p/4) there is a different cubic.
>> funs{1,5}
ans =
3.1416 3.927
>> vpa(funs{2,5},10)
ans =
0.1568856928*x^3 - 1.47861282*x^2 + 3.648107873*x - 1.731986498
With 9 data points, we will have 8 cubic segments. So the final segment has this form:
>> funs{1,8}
ans =
5.4978 6.2832
>> vpa(funs{2,8},10)
ans =
- 0.08493845396*x^3 + 1.736669488*x^2 - 10.70470091*x + 19.76765782
I could have extracted the polynomial segments in a variety of forms.
Fayyaz Ahmad
Fayyaz Ahmad 2021-10-17
clc;
clear all;
close all;
format short g
x = linspace(0,2*pi,9);
y = sin(x);
spl = spline(x,y);
bk = spl.breaks;
cf = spl.coefs;
[m,n] = size(cf);
syms z;
v = z.^((n-1:-1:0)');
eq = cf*v
------------------------------------------------------
output
----------------------------------------------------------
eq =
(4770321904107807*z)/4503599627370496 - (610766247492719*z^2)/4503599627370496 - (6120460633987311*z^3)/72057594037927936
(6206087117424961*z)/9007199254740992 + 2^(1/2)/2 - (6048313895780885*z^2)/18014398509481984 - (6120460633987189*z^3)/72057594037927936
(635448956098857*z^3)/9007199254740992 - (2413390700397705*z^2)/4503599627370496 + (159898229681981*z)/36028797018963968 + 1
(1413100694996111*z^3)/9007199254740992 - (3329540071636805*z^2)/9007199254740992 - (6365985347106939*z)/9007199254740992 + 2^(1/2)/2
(2826201389992229*z^3)/18014398509481984 - (9*z^2)/9007199254740992 - (4490500002164345*z)/4503599627370496 + 4967757600021511/40564819207303340847894502572032
(1270897912197719*z^3)/18014398509481984 + (6659080143273605*z^2)/18014398509481984 - (6365985347106939*z)/9007199254740992 - 2^(1/2)/2
- (765057579248403*z^3)/9007199254740992 + (2413390700397707*z^2)/4503599627370496 + (159898229681977*z)/36028797018963968 - 1
- (6120460633987189*z^3)/72057594037927936 + (6048313895780881*z^2)/18014398509481984 + (3103043558712477*z)/4503599627370496 - 2^(1/2)/2
  2 个评论
John D'Errico
John D'Errico 2021-10-20
编辑：John D'Errico 2021-10-20
Ran in:
@Fayyaz Ahmad
Be careful, as those polynomials cannot be used to evaluate the spline directly. (TRY IT!) You have made a mistake, if you think they can.
clc;
clear all;
close all;
format short g
x = linspace(0,2*pi,9);
y = sin(x);
spl = spline(x,y);
bk = spl.breaks;
cf = spl.coefs;
[m,n] = size(cf);
syms z;
v = z.^((n-1:-1:0)');
eq = cf*v;
For example...
eq(2)
ans = 
bk
bk = 1×9
0 0.7854 1.5708 2.3562 3.1416 3.9270 4.7124 5.4978 6.2832
So, on the interval pi/4 to pi/2, you would imply that this polynomial should approximate sin(x), correct? Try it!
P2 = matlabFunction(eq(2));
fplot(@sin,[pi/4,pi/2])
hold on
fplot(P2,[pi/4,pi/2])
legend('Sin(x)','Polynomial 2')
xlabel 'X'
grid on
Gosh. It does not look correct to me. Not even close.
Your mistake lies in a misunderstanding on your part of what spline produces, and probably why there is a problem at all. But if you will post this as an answer, you need to understand why it fails to produce something reasonable.
Fayyaz Ahmad
Fayyaz Ahmad 2021-10-20
@ John D'Errico thanks for comments, you are right. Actually the coefficients of polynomial generated by spline command are on the interval [0,bk(i+1)-bk(i)]. We have corrected the code.
-----------------------
clc ;
clear all ;
close all ;
format short g
x = linspace(pi,2*pi,9);
y = sin(x);
spl = spline(x,y);
bk = spl.breaks ;
cf = spl.coefs ;
[m,n] = size(cf );
syms z ;
v = z.^((n-1:-1:0 )');
eq = cf*v;
for i=1:m
eq(i) = simplify(subs(eq(i),z,z-bk(i)));
end
figure(1)
hold on
q = length(bk);
for i=1:q-1
s1 = linspace(bk(i),bk(i+1),10);
p = matlabFunction(eq(i));
fplot(@sin,[bk(i),bk(i+1)])
plot(s1,p(s1),'.-r')
pause
end
eq
--------------------------------------
output is
--------------------------------------
eq =
(2260863087144877*pi)/2251799813685248 - (2260863087144877*z)/2251799813685248 + (661015214105797*(z - pi)^2)/36028797018963968 + (651968472240021*(z - pi)^3)/4503599627370496 + 4967757600021511/40564819207303340847894502572032
(9349213407344919*pi)/9007199254740992 - (1038801489704991*z)/1125899906842624 + (1701418325597523*(z - (9*pi)/8)^2)/9007199254740992 + (2607873888959953*(z - (9*pi)/8)^3)/18014398509481984 - 6893811853601121/18014398509481984
(15927176730973595*pi)/18014398509481984 - (3185435346194719*z)/4503599627370496 - 2^(1/2)/2 + (6475165695337053*(z - (5*pi)/4)^2)/18014398509481984 + (6610892248307199*(z - (5*pi)/4)^3)/72057594037927936
(9475875933393265*pi)/18014398509481984 - (861443266672115*z)/2251799813685248 + (4211117090838325*(z - (11*pi)/8)^2)/9007199254740992 + (4785409777295905*(z - (11*pi)/8)^3)/144115188075855872 - 2080391759176529/2251799813685248
(27*pi)/144115188075855872 - (9*z)/72057594037927936 + (570433996317983*(z - (3*pi)/2)^2)/1125899906842624 - (2392704888647973*(z - (3*pi)/2)^3)/72057594037927936 - 1
(3445773066688457*z)/9007199254740992 - (44795049866949941*pi)/72057594037927936 + (8422234181676639*(z - (13*pi)/8)^2)/18014398509481984 - (6610892248307199*(z - (13*pi)/8)^3)/72057594037927936 - 8321567036706117/9007199254740992
(6370870692389435*z)/9007199254740992 - (44596094846726045*pi)/36028797018963968 - 2^(1/2)/2 + (809395711917129*(z - (7*pi)/4)^2)/2251799813685248 - (5215747777920001*(z - (7*pi)/4)^3)/36028797018963968
(4155205958819961*z)/4503599627370496 - (62328089382299415*pi)/36028797018963968 + (6805673302390011*(z - (15*pi)/8)^2)/36028797018963968 - (5215747777919999*(z - (15*pi)/8)^3)/36028797018963968 - 3446905926800567/9007199254740992

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Spline Postprocessing 的更多信息

标签

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by