# -*- coding: cp949 -*-
def __newton_subcalc(sx, sy, ex, ey):
return (ey - sy)/float(ex - sx) # 기울기를 구한다
def __newton_build(srcList, stack):
"""
x y
1 11
2 12 F[1,2] = (1, 11, 2, 12)
3 13 F[2,3] = (2, 12, 3, 13) F[1,3] = (1, F[1,1], 3, F[2,3])
"""
stack.append((srcList[0][2], srcList[-1][2]))
if len(srcList) <= 1:
return stack
sx, mx, sy = srcList[0]
dstList = []
for mx, ex, ey in srcList[1:]:
dstList.append((sx, ex, __newton_subcalc(sx, sy, ex, ey)))
sx, sy = mx, ey
return __newton_build(dstList, stack)
def newton_interpolation(srcList, x):
"""
stack
P(x) = stack[0] + stack[1](x - x0)
+ stack[2](x - x0)(x - x1)
+ stack[3](x - x0)(x - x1)(x - x2)
...
x 가 리스트 앞쪽에 가까울 때는 전향차분을 사용하고 (왼쪽 리턴값)
x 가 리스트 뒤쪽에 가까울 때는 후향차분을 사용한다. (오른쪽 리턴값)
"""
# 재귀함수를 호출할때 사용하는 형식으로 데이터를 바꾼다
srcList = [(cx, cx, cy) for cx, cy in srcList]
stack = __newton_build(srcList, [])
fc, bc = stack[0]
fr = fc
count = 1
for fc, bc in stack[1:]:
cx, ct, cy = srcList[0]
dx = (x - cx)
for cx, ct, cy in srcList[1:count]:
dx *= (x - cx)
fr += fc * dx
count += 1
revList = [] + srcList
revList.reverse()
fc, bc = stack[0]
br = bc
count = 1
for fc, bc in stack[1:]:
cx, ct, cy = revList[0]
dx = (x - cx)
for cx, ct, cy in revList[1:count]:
dx *= (x - cx)
br += bc * dx
count += 1
return fr, br
positions = [
(1.0, 0.7651977),
(1.3, 0.6200860),
(1.6, 0.4554022),
(1.9, 0.2818186),
(2.2, 0.1103623),
]
print newton_interpolation(positions, 1.5)

(
0)

(
0)
http://imp17.com/tc/myevan/trackback/232