Script to solve second order (polynomial) interpolation | Unix Linux Forums | Shell Programming and Scripting

 Member Name Remember Me? Password

 Shell Programming and Scripting Post questions about KSH, CSH, SH, BASH, PERL, PHP, SED, AWK and OTHER shell scripts and shell scripting languages here.

# Script to solve second order (polynomial) interpolation

## Shell Programming and Scripting

 Tags awk, numbering, text processing

#1
01-21-2013
 Tzeronone Registered User Join Date: Jan 2013 Last Activity: 28 July 2013, 10:06 PM EDT Posts: 10 Thanks: 6 Thanked 0 Times in 0 Posts
Script to solve second order (polynomial) interpolation

Currently I have awk command to do linear interpolation

Code:
```awk '
{
P[\$1]=\$2
I[i++]=\$1
}
END {
j=0; s=I[j]; t=I[j+1]
for(i=m;i<=n;i++) {
if(I[j+2] && i>t) {
j++; s=I[j]; t=I[j+1]
}
print i, P[s]+(i-s)*(P[t]-P[s])/(t-s)
}
} ' m=1 n=8 infile```

FILE CONTENT INPUT# a.txt

Code:
```#X Y
1 22.3125
4 22.5
8 22.1875```

Any idea how to change the code to polynomial interpolation?
Assume that Y of X0 = 0. Please help me...

Last edited by Scrutinizer; 01-21-2013 at 03:54 AM.. Reason: icode tags => code tags
Sponsored Links
#2
01-21-2013
 Corona688   Mead Rotor Join Date: Aug 2005 Last Activity: 17 April 2014, 6:29 PM EDT Location: Saskatchewan Posts: 18,531 Thanks: 681 Thanked 3,036 Times in 2,856 Posts
What do you mean by polynomial interpolation?
Sponsored Links
#3
01-21-2013
 Tzeronone Registered User Join Date: Jan 2013 Last Activity: 28 July 2013, 10:06 PM EDT Posts: 10 Thanks: 6 Thanked 0 Times in 0 Posts
Quote:
 Originally Posted by Corona688 What do you mean by polynomial interpolation?

Polynomial interpolation is an interpolation that based on three value points (two previous points and next point). So, if we see the data, there are an ID 1, 4, and 8. So, we need to find out the value of 2, 3, 5, 6, and 7 based on value 1, 4, and 8. The formula is that
(((x-x2) * (x-x3)) / ((x1-x2) * (x1-x3))) * y1 + (((x-x1) * (x-x3)) / ((x2-x1) * (x2-x3))) * y2 + (((x-x1) * (x-x2)) / ((x3-x1) * (x3-x1))) * y3
x = current ID;
x1 = the first known ID (second previous known ID); ---> 1
x2 = the second known ID (first previous known ID); ---> 4
x3 = the third known ID (next known ID); ---> 8
y1 = the first known value (the value of ID x1)
y2 = the second known value (the value of ID x2)
y3 = the third known value (the value of ID x3)
#4
01-21-2013
 Chubler_XL   Registered User Join Date: Oct 2010 Last Activity: 17 April 2014, 12:50 PM EDT Posts: 2,601 Thanks: 94 Thanked 804 Times in 757 Posts
How about this:

Code:
```awk '
{ P[\$1]=\$2 ; m=\$1>m?\$1:m; x[NR]=\$1 ; y[NR]=\$2 }
END {
for(i=1;i<=m;i++) {
if (i in P) printf "%0.4f %0.4f\n", i, P[i];
else printf "%0.4f %0.4f\n", i, \
(((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1] + \
(((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2] + \
(((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[1]))) * y[3] ;
}
}' infile```

 The Following User Says Thank You to Chubler_XL For This Useful Post: Tzeronone (01-21-2013)
Sponsored Links
#5
01-21-2013
 Tzeronone Registered User Join Date: Jan 2013 Last Activity: 28 July 2013, 10:06 PM EDT Posts: 10 Thanks: 6 Thanked 0 Times in 0 Posts
Quote:
 Originally Posted by Chubler_XL How about this: Code: ```awk ' { P[\$1]=\$2 ; m=\$1>m?\$1:m; x[NR]=\$1 ; y[NR]=\$2 } END { for(i=1;i<=m;i++) { if (i in P) printf "%0.4f %0.4f\n", i, P[i]; else printf "%0.4f %0.4f\n", i, \ (((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1] + \ (((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2] + \ (((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[1]))) * y[3] ; } }' infile```
It is working but shows the wrong result, thus it remains dissolved. Any how, thanks.
Sponsored Links
#6
01-21-2013
 Chubler_XL   Registered User Join Date: Oct 2010 Last Activity: 17 April 2014, 12:50 PM EDT Posts: 2,601 Thanks: 94 Thanked 804 Times in 757 Posts
Damn,

Here is a debug version that show the formula with expanded values, it might help us find the issue.

Note data file should not contain the heading ( #X Y ) line:

Code:
```awk '
{ P[\$1]=\$2 ; m=\$1>m?\$1:m; x[NR]=\$1 ; y[NR]=\$2 }
END {
for(i=1;i<=m;i++) {
if (i in P) printf "%0.4f %0.4f\n", i, P[i];
else { printf "%0.4f %0.4f\n", i, \
(((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1] + \
(((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2] + \
(((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[1]))) * y[3] ;
printf "   (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f +\n",
i,x[2],i,x[3],x[1],x[2],x[1],x[3],y[1];
printf "   (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f +\n",
i,x[1],i,x[3],x[2],x[1],x[2],x[3],y[2];
printf "   (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f\n",
i,x[1],i,x[2],x[3],x[1],x[3],x[1],y[3];
printf "  = %f + %f + %f\n",
(((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1],
(((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2],
(((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[1]))) * y[3] ;
}
}
}' infile```

Exampe output for X=7.0:

Code:
```7.0000 16.2130
(((7.000000-4.000000) * (7.000000-8.000000)) / ((1.000000-4.000000) * (1.000000-8.000000))) * 22.312500 +
(((7.000000-1.000000) * (7.000000-8.000000)) / ((4.000000-1.000000) * (4.000000-8.000000))) * 22.500000 +
(((7.000000-1.000000) * (7.000000-4.000000)) / ((8.000000-1.000000) * (8.000000-1.000000))) * 22.187500
= -3.187500 + 11.250000 + 8.150510```

Last edited by Chubler_XL; 01-21-2013 at 11:24 PM.. Reason: Fix indenting
 The Following User Says Thank You to Chubler_XL For This Useful Post: Tzeronone (01-22-2013)
Sponsored Links
#7
01-22-2013
 Tzeronone Registered User Join Date: Jan 2013 Last Activity: 28 July 2013, 10:06 PM EDT Posts: 10 Thanks: 6 Thanked 0 Times in 0 Posts
Quote:
 Originally Posted by Chubler_XL Damn, Here is a debug version that show the formula with expanded values, it might help us find the issue. Note data file should not contain the heading ( #X Y ) line: Code: ```awk ' { P[\$1]=\$2 ; m=\$1>m?\$1:m; x[NR]=\$1 ; y[NR]=\$2 } END { for(i=1;i<=m;i++) { if (i in P) printf "%0.4f %0.4f\n", i, P[i]; else { printf "%0.4f %0.4f\n", i, \ (((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1] + \ (((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2] + \ (((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[1]))) * y[3] ; printf " (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f +\n", i,x[2],i,x[3],x[1],x[2],x[1],x[3],y[1]; printf " (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f +\n", i,x[1],i,x[3],x[2],x[1],x[2],x[3],y[2]; printf " (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f\n", i,x[1],i,x[2],x[3],x[1],x[3],x[1],y[3]; printf " = %f + %f + %f\n", (((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1], (((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2], (((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[1]))) * y[3] ; } } }' infile``` Exampe output for X=7.0: Code: ```7.0000 16.2130 (((7.000000-4.000000) * (7.000000-8.000000)) / ((1.000000-4.000000) * (1.000000-8.000000))) * 22.312500 + (((7.000000-1.000000) * (7.000000-8.000000)) / ((4.000000-1.000000) * (4.000000-8.000000))) * 22.500000 + (((7.000000-1.000000) * (7.000000-4.000000)) / ((8.000000-1.000000) * (8.000000-1.000000))) * 22.187500 = -3.187500 + 11.250000 + 8.150510```
GREAT. It is working. But how to declare if x=0, then y[x]=0. Because for y[3] and y[4], the x[1]=0, thus y[0]=0.

Code:
```awk '
{ P[\$1]=\$2 ; m=\$1>m?\$1:m; x[NR]=\$1 ; y[NR]=\$2 }
END {
for(i=1;i<=m;i++) {
if (i in P) printf "%0.4f %0.4f\n", i, P[i];
else { printf "%0.4f %0.4f\n", i, \
(((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1] + \
(((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2] + \
(((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[2]))) * y[3] ;
printf "   (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f +\n",
i,x[2],i,x[3],x[1],x[2],x[1],x[3],y[1];
printf "   (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f +\n",
i,x[1],i,x[3],x[2],x[1],x[2],x[3],y[2];
printf "   (((%f-%f) * (%f-%f)) / ((%f-%f) * (%f-%f))) * %f\n",
i,x[1],i,x[2],x[3],x[1],x[3],x[2],y[3];
printf "  = %f + %f + %f\n",
(((i-x[2]) * (i-x[3])) / ((x[1]-x[2]) * (x[1]-x[3]))) * y[1],
(((i-x[1]) * (i-x[3])) / ((x[2]-x[1]) * (x[2]-x[3]))) * y[2],
(((i-x[1]) * (i-x[2])) / ((x[3]-x[1]) * (x[3]-x[2]))) * y[3] ;
}
}
}' infile```

Last edited by Tzeronone; 01-22-2013 at 02:00 AM.. Reason: icode tags changed to code tags
Sponsored Links

 Thread Tools Search this Thread Search this Thread: Advanced Search Display Modes Linear Mode

 More UNIX and Linux Forum Topics You Might Find Helpful Thread Thread Starter Forum Replies Last Post nex_asp Shell Programming and Scripting 12 12-05-2012 08:38 AM ardy_yana Shell Programming and Scripting 7 04-05-2012 03:04 AM z1dane Web Programming 0 01-27-2010 12:10 AM egkumpe UNIX for Dummies Questions & Answers 0 08-05-2004 06:05 PM

All times are GMT -4. The time now is 05:32 PM.