Commit e4813518 authored by Thomas Bruns's avatar Thomas Bruns
Browse files

four parameter fit added (single and sequentiel)

parent f528a158
......@@ -168,8 +168,8 @@ def seq_threeparsinefit(y,t,f0):
"""
Tau = 1.0/f0
dt = t[1]-t[0]
N = Tau/dt ## samples per section
M = sp.floor(t.size/N) ## number of sections or periods
N = int(Tau/dt) ## samples per section
M = int(sp.floor(t.size/N)) ## number of sections or periods
abc = sp.zeros((M,3))
......@@ -179,9 +179,88 @@ def seq_threeparsinefit(y,t,f0):
abc[i,:] = (threeparsinefit(yi,ti,f0))
return abc ## matrix of all fit vectors per period
#def fourparsinefit(y,t,f0,df_max=None):
# if dy_max is None :
# df_max = 0.05*f0
# four parameter sine-fit (with frequency approximation)
def fourparsinefit(y,t,f0, tol=1.0e-7, nmax=1000):
"""
y sampled data values \n
t sample times of y \n
f0 estimate of sine frequency \n
tol rel. frequency correction where we stop \n
nmax maximum number of iterations taken \n
\n
returns the vector [a, b, c, w] of a*sin(w*t)+b*cos(w*t)+c
"""
abcd = threeparsinefit(y,t,f0)
w = 2*sp.pi*f0
err = 1
i = 0
while ((err > tol) and (i<nmax)):
D = sp.array([sp.sin(w*t),
sp.cos(w*t),
sp.ones(t.size),
abcd[0]*t*sp.cos(w*t)-abcd[1]*t*sp.sin(w*t)
])
abcd = (la.lstsq(D.transpose(),y))[0]
dw = abcd[3]
w = w+0.9*dw
err = dw/w
assert i < nmax, "iteration error"
return sp.hstack((abcd[0:3],w/(2*sp.pi)))
"""
from octave ...
function abcw = fourParSinefit(data,w0)
abc = threeParSinefit(data,w0);
a=abc(1);
b=abc(2);
c=abc(3);
w = w0;
do
D = [sin(w.*data(:,1)) , cos(w.*data(:,1)) , ones(rows(data),1) , a.*data(:,1).*cos(w.*data(:,1)) - b.*data(:,1).*sin(w.*data(:,1)) ];
s = D \ data(:,2);
dw = s(4);
w = w+0.9*dw;
err = abs(dw/w);
until (err < 1.0e-8 );
abcw = [s(1),s(2),s(3),w];
endfunction
"""
# periodical sinefit at known frequency
def seq_fourparsinefit(y,t,f0, tol=1.0e-7, nmax=1000):
"""
period-wise sine-fit at a known frequency\n
y vector of sample values \n
t vector of sample times\n
f0 known frequency\n
\n
returns a (n,3)-matrix of coefficient-triplets [[a,b,c], ...]\n
for y = a*sin(2*pi*f0*t) + b*cos(2*pi*f0*t) + c
"""
Tau = 1.0/f0
dt = t[1]-t[0]
N = int(sp.floor(Tau/dt)) ## samples per section
M = int(sp.floor(t.size/N)) ## number of sections or periods
abcd = sp.zeros((M,4))
for i in range(M):
ti = t[i*N:(i+1)*N]
yi = y[i*N:(i+1)*N]
abcd[i,:] = (fourparsinefit(yi,ti,f0,tol=tol,nmax=nmax))
return abcd ## matrix of all fit vectors per period
# fitting a pseudo-random multi-sine signal with 2*Nf+1 parameters
def multi_threeparsinefit(y,t,f0): # fo vector of frequencies
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment