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

general check and clean up

parent ddd35192
......@@ -428,7 +428,7 @@ def seq_threeparcounterfit(y, t, f0, diff=False):
if diff=True use differentiation to remove carrier (c.f. source)
"""
Tau = 1.0 / f0
dt = t[1] - t[0]
dt = np.mean(np.diff(t))
N = int(np.floor(Tau / dt)) ## samples per section
M = int(np.floor(t.size / N)) ## number of sections or periods
......@@ -731,7 +731,7 @@ def seq_multi_threeparam_Dmatrix(f,t,periods=1, progressive=True):
ci = []
ri = []
fr = []
# Designmatrix for sin/cos
# Designmatrix for cos/sin
for fi, omi in zip(f,2*np.pi*f):
Nri = 0 # counter for the current row index
if progressive:
......@@ -740,15 +740,16 @@ def seq_multi_threeparam_Dmatrix(f,t,periods=1, progressive=True):
tau = 1/fi*periods # slice length in seconds for periods of fi
t_sl = np.array_split(t,np.ceil(T/tau)) # array of slices of sample times
fr = fr + [fi]*len(t_sl) # len(t_sl) times frequency fi in Design matrix
for ti in t_sl: # loop over slices fo one frequncy
sin = np.sin(omi * ti)
for ti in t_sl: # loop over slices for one frequency
# cosine part
cos = np.cos(omi * ti)
data = np.hstack((data,cos)) # data vector
ci = ci+[col]*len(ti) # column index vector
col = col+1
ri = ri+ [Nri+i for i in range(len(ti))] # row index vector
#sine part
sin = np.sin(omi * ti)
data = np.hstack((data,sin)) # data vector
ci = ci+[col]*len(ti) # column index vector
col = col+1
......@@ -756,7 +757,7 @@ def seq_multi_threeparam_Dmatrix(f,t,periods=1, progressive=True):
Nri = Nri+len(ti)
# Designmatrix for Bias
# Bias part
data = np.hstack((data,np.ones((len(t)))))
ci = ci = ci+[col]*len(t)
ri = ri + [i for i in range(len(t))]
......@@ -794,8 +795,8 @@ def seq_multi_threeparsinefit(f,y,t,periods=1, D_fr=None, abc0=None, progressive
if D_fr is None:
fr, D = seq_multi_threeparam_Dmatrix(f,t,periods, progressive=progressive) # calculate the design matrix (as sparse matrix)
else:
D = D_fr[0]
fr = D_fr[1]
D = D_fr[0] # vector of sparse design matrix
fr = D_fr[1] # vector of frequencies
abc = sp_lsqr(D,y,x0=abc0, atol=1.0e-9, btol=1.0e-9)
y0 = D.dot(abc[0])
......@@ -803,10 +804,10 @@ def seq_multi_threeparsinefit(f,y,t,periods=1, D_fr=None, abc0=None, progressive
#print(abc)
f_ab_c =[]
k=0
for fi in fr:
f_ab_c = f_ab_c+ [fi, abc[0][k],abc[0][k+1]]
for fi in fr: # compile a list of frequencies and coefficients
f_ab_c = f_ab_c+ [fi, abc[0][k],abc[0][k+1]] # [f, a, b]...
k=k+2
f_ab_c = f_ab_c + [0.0,abc[0][-1],0.0]
f_ab_c = f_ab_c + [0.0,abc[0][-1],0.0] # add the bias to the list [0,c,0]
f_ab_c = np.array(f_ab_c).reshape((len(f_ab_c)//3,3))
#print(f_ab_c)
......@@ -895,8 +896,8 @@ def seq_multi_amplitude(f_ab_c):
"""
#print("Test")
N = f_ab_c.shape[0]-1
amp = np.array([1j*np.abs(f_ab_c[i,1]+f_ab_c[i][2]) for i in range(f_ab_c.shape[0]-1)]).reshape((N,1))
N = f_ab_c.shape[0]-1 # exclude the bias
amp = np.array([np.abs(1j*f_ab_c[i,1]+f_ab_c[i][2]) for i in range(f_ab_c.shape[0]-1)]).reshape((N,1))
return np.hstack((f_ab_c[:-1,0].reshape((N,1)), amp))
def seq_multi_phase(f_ab_c, deg=True):
......@@ -907,8 +908,9 @@ def seq_multi_phase(f_ab_c, deg=True):
Parameters
----------
f_ab_c : 2-d numpy array of floats (Nx3)
f_ab_c : 2-d numpy array of floats (Nx3)
f,a,b in a row for several rows, as returned by seq_multi_threeparsinefit.
x=a*cos+b*sin+c
deg : Boolean, optional
Flag whether result is in degrees or rad. The default is True (Degrees).
......@@ -918,8 +920,8 @@ def seq_multi_phase(f_ab_c, deg=True):
frequency and associated initial phase.
"""
N = f_ab_c.shape[0]-1
phase = np.array([1j*np.angle(f_ab_c[i,1]+f_ab_c[i][2],deg=deg) for i in range(f_ab_c.shape[0]-1)]).reshape((N,1))
N = f_ab_c.shape[0]-1 # exclude the bias
phase = np.array([np.angle(1j*f_ab_c[i,1]+f_ab_c[i][2],deg=deg) for i in range(f_ab_c.shape[0]-1)]).reshape((N,1))
return np.hstack((f_ab_c[:-1,0].reshape((N,1)), phase))
def seq_multi_bias(f_ab_c):
......
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