Re: [討論] fit平面問題

看板MATLAB作者tn00364361 (小氵斿)時間10年前 (2015/08/18 20:07)推噓1(1推 0噓 6→)

留言7則, 2人參與討論串2/2 (看更多)

※ 引述《liwes5566 (wes5566)》之銘言： : 目前是利用指令中的fit(2元2次多項式poly22)下去fit : 再一般的情況下，都可以fit得很好 : 可是當我是fit一個平面(數值全部都一樣)，所得到的答案不是平面 : 這是我已1024*1024矩陣，數值為2^16的結果 : fitobject(x,y)=p00+p10x+p01y+p20x^2+p11xy+p02y^2 : p00=6.554e+04 : p10=-5.012e-12 : p01=3.291e-09 : p20=2.868e-14 : p11=-7.413e-16 : p02=-1.801e-11 : 有沒有甚麼方法可以讓p00的值為2^16，其他都為0呢? 任何方法fit的結果都不太可能會是真的零誤差在一般電腦上MATLAB的machine epsilon大概是2.22e-16 通常可接受的數值誤差範圍大概是sqrt(machine epsilon)，也就是1.49e-8 ---------------------------------------- Fit多項式曲面(或是曲線)可以看成是一個線性的least-squares問題假設你有m個資料點，p = [p1; p2; ...; pn]為你想要fit的曲面的係數。而曲面的方程式可以寫成 z = a1 * p1 + a2 * p2 + ... + an * pn = dot(a, p) 其中a是x跟y的函數把所有的資料點"疊"起來，可以寫成 z1 ~ a11 * p1 + a12 * p2 + ... + a1n * pn z2 ~ a21 * p1 + a22 * p2 + ... + a2n * pn ... zm ~ am1 * p1 + am2 * p2 + ... + amn * pn 這可以重新寫成矩陣的形式: Z ~ A * p (這裡的"~"是指大約等於) 其中，[A]ij = aij，Z = [z1; z2; ...; zm] 這時候的最佳解(依least-squares)就是 p_fit = argmin(norm(A * p - Z)) 也就是 p_fit = A \ Z (或是寫成 p_fit = pinv(A) * Z) ---------------------------------------- 回到你的case。先把格點跟資料點都向量化，變成column vector n = 1024; x = reshape(x, n^2, 1); % grid point y = reshape(y, n^2, 1); % grid point z = reshape(z, n^2, 1); % data point ---------------------------------------- 1. 如果用least-squares去fit所有的係數: A = [ones(n^2, 1), x, y, x.^2, x .* y, y.^2]; p_fit = A \ z; % = [p00; p10; p01; p20; p11; p02] ---------------------------------------- 2. 如果只去fit除了p00以外的係數(已知p00 = 2^16): z = z - 2^16; A = [x, y, x.^2, x .* y, y.^2]; p_fit = A \ z; % = [p10; p01; p20; p11; p02] ---------------------------------------- 這個方法可以很輕易的fit任意維度空間的任意多項式曲面(例如在三維空間中fit橢球) 希望這對你有幫助 -- ※ 發信站: 批踢踢實業坊(ptt.cc), 來自: 68.43.178.167 ※ 文章網址: https://www.ptt.cc/bbs/MATLAB/M.1439899632.A.F90.html ※ 編輯: tn00364361 (68.43.178.167), 08/18/2015 22:09:36