%求取導數公式的系數
m=input('m=');					%鍵盤輸入m
k=input('k=');					%鍵盤輸入k，k-1次多項式
p=input('p=');                  %鍵盤輸入實驗數據點數p
n=2*m+1;
y=[0.142,0.176,0.215,0.260,0.309,0.363,0.419,0.479,0.539,0.598,0.654,0.706,0.752,0.789,0.817,0.834,0.840,0.834,0.817,0.789,0.752,0.706,0.654,0.598,0.539,0.479,0.419,0.363,0.309,0.260,0.215,0.176,0.142,0.113,0.089,0.068,0.052,0.039,0.029,0.021];
x=zeros(n,k);
for i=-m:m							%構造X矩陣行循環
    for j=0:k-1						%構造X矩陣列循環
	    x(i+m+1,j+1)=i^j;		
    end
end
d=inv(x'*x)*x';                     %導數公式系數
yc=zeros(k - 1,p - 2 * m);          %存放導數結果的矩陣
ax=1:1:p;                           %原始數據點的序號
axc=m+1:1:p-m;                      %導數結果的序號
for i = 2:k
    for j = m + 1:p - m
        for q = 1: n
            yc(i-1,j-m)=yc(i-1,j-m)+y(j-m+q-1)*d(i,q);      %導數計算
        end
    end
end
subplot(1,2,1)                      %設置畫圖區域為水平排列,兩個區域
plot(ax,y,'.-',axc,yc,'.-r')        %畫出原始數據與導數結果的圖形
subplot(1,2,2)                      %設置圖形區域
plot(axc,yc,'.-r')                  %畫出導數結果圖