lagrange interpolation, .m

can anyone explain me how to use this program
function y=lagrange(x,pointx,pointy)
%
%LAGRANGE approx a point-defined function using the Lagrange polynomial interpolation
%
% LAGRANGE(X,POINTX,POINTY) approx the function definited by the points:
% P1=(POINTX(1),POINTY(1)), P2=(POINTX(2),POINTY(2)), ..., PN(POINTX(N),POINTY(N))
% and calculate it in each elements of X
%
% If POINTX and POINTY have different number of elements the function will return the NaN value
%
% function wrote by: Calzino
% 7-oct-2001
%
n=size(pointx,2);
L=ones(n,size(x,2));
if (size(pointx,2)~=size(pointy,2))
fprintf(1,'\nERROR!\nPOINTX and POINTY must have the same number of elements\n');
y=NaN;
else
for i=1:n
for j=1:n
if (i~=j)
L(i,:)=L(i,:).*(x-pointx(j))/(pointx(i)-pointx(j));
end
end
end
y=0;
for i=1:n
y=y+pointy(i)*L(i,:);
end
end

2 Comments

How to run this code
Hardipsinh Jadeja
Hardipsinh Jadeja on 24 Apr 2018
Edited: Hardipsinh Jadeja on 24 Apr 2018
If size of pointx and pointy is same size then why not print the statement

Sign in to comment.

 Accepted Answer

pointx and pointy are two vectors of data values, x is a vector of points where you want to interpolate. For example:
x = 0:10;
y = x.^2;
xx = linspace(0,10);
yy = lagrange(xx,x,y);
plot(x,y,'o',xx,yy,'.')
As an aside, with no offense intended to Calzino, there are other options available for interpolation. Firstly, of course, interp1 is a standard MATLAB function, with options for linear, cubic spline, and PCHIP interpolation. Cleve Moler (aka The Guy Who Wrote MATLAB) also has a Lagrange interpolation function available for download.

7 Comments

buxZED
buxZED on 17 Mar 2011
still dont get it
given the points x0 , x1 , x2 , x3 , x4
and the values f(x0 ), f(x1 ), f(x2 ), f(x3 ), f(x4 )
how can i get an 4th order accurat function f in position x
just point me in the right derection
Jan
Jan on 17 Mar 2011
I think, Matt Tearle's example does point in the right direction already. But I can repeat it:
y = lagrange(x, [x0,x1,x2,x3,x4], [f(x0),f(x1),f(x2),f(x3),f(x4)]);
This is exactly found in the help section of the function.
Right, what Jan said. In my example,x and y are vectors of the points x0, x1, ..., x4 and f(x0), ..., f(x4). The new point you're calling x is what I called xx. I used a vector of points, but it could be a single value.
"Jan Simon on 17 Mar 2011 I think, Matt Tearle's example does point in the right direction already. But I can repeat it: y = lagrange(x, [x0,x1,x2,x3,x4], [f(x0),f(x1),f(x2),f(x3),f(x4)]); This is exactly found in the help section of the function." I gave some values like , y = lagrange(x, [1,2,3,4], [2,4,6,8]) which returned, 2 4 6 8 10
what does this means?
Jan
Jan on 12 Apr 2017
@Mudra: What does what mean? Please open a new thread for a new question. Then provide ay many details as required to repdoduce the problem.
@Matt Tearle
That link is dead, I don't suppose you have an updated one?
Note: the File Exchange has some more advanced polyinterp functions.

Sign in to comment.

More Answers (5)

Matt Fig
Matt Fig on 16 Mar 2011

1 vote

This is really a question for the author of the program. I believe it is also bad etiquette to post somebody's code like that without permission.
Did you try to contact the author?

3 Comments

It's from File Exchange, so I don't seem any great harm in posting it.
Matt Fig
Matt Fig on 16 Mar 2011
Ah, but I wasn't talking about harm, just polite behavior. The author should have been contacted first, that's all.
Fair call. I guess it does open the door for people to bash the author's code in a separate location, which would be uncool.

Sign in to comment.

%% Lagrangian interpolation
clear;clc;close all;
X=[-3 -2.5 -1 0 2 3.75 4.25 7];
Y=(sqrt(1+abs(X)));
xq=min(X):0.1:max(X);
f=(sqrt(1+abs(xq)));
syms x
S=0;
for i=1:length(X)
temp=X;
A=temp(i);
temp(i)=[];
L=prod((x-temp)./(A-temp),'all');
S=(L*Y(i))+S;
L=[];
end
figure()
fplot(S,'black--',[min(X) max(X)]);
hold on
F=interp1(X,Y,xq);
plot(xq,F,"bo");
hold on
plot(xq,f,"r*");
legend("Lagrangian","interp1","f(x)",'Location','north');
xlabel(" X axis ");
ylabel(" Y axis");
title("Lagrangian interpolation VS interp1-MatlabFunction")
Above we can see an easy way to implement lagrangian interpolation which has been checked with matlab interp1() function;
From MohammadReza Arani
mohammadrezaarani@ut.ac.ir

4 Comments

how to imaplement this method for image?
X = 1:numel(YourImage);
Y = double(YourImage(:));
This is not likely to give good results.
Perhaps you want something different than this? Perhaps the image has a curve already drawn in it, and you want to extract the points from the image and then do interpolation on the curve?
yes, but how to code all this?
See https://www.mathworks.com/matlabcentral/fileexchange/?term=tag:%22digitize%22 for a number of File Exchange contributions that try to extract data from images of graphs.

Sign in to comment.

norah
norah on 10 May 2023

0 votes

how can i find error bound ?
John
John on 31 Jul 2023
Edited: Walter Roberson on 10 Oct 2023
function Y = Lagrange_371(x,y,X)
n = length(x) - 1;
Y = 0;
for i = 0:n
prod = 1;
for j = 0:n
if i ~= j
prod = prod.*(X - x(j+1))./(x(i+1) - x(j+1));
end
end
Y = Y + prod*y(i+1);
end
end

1 Comment

bonsoir ,comment appliquer cette fonction sur un exemple

Sign in to comment.

MUHAMMAD IQRAM HAFIZ
MUHAMMAD IQRAM HAFIZ on 21 May 2024

0 votes

function P = lagrange_interpolation_3point(x1, y1, x2, y2, x3, y3, x)
% Compute the Lagrange basis polynomials
L1 = ((x - x2) .* (x - x3)) / ((x1 - x2) * (x1 - x3));
L2 = ((x - x1) .* (x - x3)) / ((x2 - x1) * (x2 - x3));
L3 = ((x - x1) .* (x - x2)) / ((x3 - x1) * (x3 - x2));
% Compute the Lagrange polynomial
P = y1 * L1 + y2 * L2 + y3 * L3;
end
% Given points
x0 = 0; y0 = 0;
x1 = 2; y1 = 10;
x2 = 4; y2 = 20;
x3 = 8; y3 = 100;
% Point at which to evaluate the polynomial
x = 4.5;
% Calculate the interpolation polynomial at x = 4.5
P = lagrange_interpolation_3point(x1, y1, x2, y2, x3, y3, x);
% Display the result
fprintf('The interpolated value at x = %.1f is P = %.2f\n', x, P);

Categories

Community Treasure Hunt

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

Start Hunting!