|
@@ -1,4 +1,4 @@
|
|
-from cmath import atan
|
|
|
|
|
|
+# from cmath import atan
|
|
import numpy as np
|
|
import numpy as np
|
|
|
|
|
|
|
|
|
|
@@ -14,12 +14,13 @@ def open_file(path):
|
|
re.append(float(temp[1]))
|
|
re.append(float(temp[1]))
|
|
im.append(float(temp[2]))
|
|
im.append(float(temp[2]))
|
|
return freq, re, im
|
|
return freq, re, im
|
|
-
|
|
|
|
|
|
+
|
|
|
|
+
|
|
def prepare_data(freq, re, im, fl=None):
|
|
def prepare_data(freq, re, im, fl=None):
|
|
"""the function takes raw data and gives vectors of eq (8)"""
|
|
"""the function takes raw data and gives vectors of eq (8)"""
|
|
# finding fl from the point with smallest magnitude if argument not provided
|
|
# finding fl from the point with smallest magnitude if argument not provided
|
|
if fl is None:
|
|
if fl is None:
|
|
- s = abs(np.array(re) + np.array(im)*1j)
|
|
|
|
|
|
+ s = abs(np.array(re) + np.array(im) * 1j)
|
|
# frequency of loaded resonance
|
|
# frequency of loaded resonance
|
|
fl = freq[list(abs(s)).index(min(abs(s)))]
|
|
fl = freq[list(abs(s)).index(min(abs(s)))]
|
|
|
|
|
|
@@ -39,6 +40,7 @@ def prepare_data(freq, re, im, fl=None):
|
|
data = np.array([e1, e2, e3, gamma, p], dtype=np.cdouble)
|
|
data = np.array([e1, e2, e3, gamma, p], dtype=np.cdouble)
|
|
return data, fl
|
|
return data, fl
|
|
|
|
|
|
|
|
+
|
|
def solution(data):
|
|
def solution(data):
|
|
""" takes projections of equation (8) on vectors e1, e2, e3 and solves the equations.
|
|
""" takes projections of equation (8) on vectors e1, e2, e3 and solves the equations.
|
|
It is also gives matrix of equation"""
|
|
It is also gives matrix of equation"""
|
|
@@ -68,27 +70,39 @@ def q_factor(a):
|
|
def recalculation_of_data(data, a, c, d, error=False):
|
|
def recalculation_of_data(data, a, c, d, error=False):
|
|
"""preparation data for the next iteration of solving system"""
|
|
"""preparation data for the next iteration of solving system"""
|
|
# data = np.array([e1, e2, e3, gamma, p], dtype=complex), t = e1, 1 = e2
|
|
# data = np.array([e1, e2, e3, gamma, p], dtype=complex), t = e1, 1 = e2
|
|
- eps = np.array(a[0]*data[0] + a[1]*data[1] - a[2]*data[0]*data[3] - data[3], dtype=complex)
|
|
|
|
|
|
+ eps = np.array(a[0] * data[0] + a[1] * data[1] - a[2] * data[0] * data[3] - data[3], dtype=complex)
|
|
# eps is eq(7) line's errors
|
|
# eps is eq(7) line's errors
|
|
- S2 = np.dot(abs(data[4]), abs(eps)*abs(eps)) # the weighted squared sum of errors
|
|
|
|
|
|
+ S2 = np.dot(abs(data[4]), abs(eps) * abs(eps)) # the weighted squared sum of errors
|
|
sigma2A = [] # the square of standart deviation coefficients a
|
|
sigma2A = [] # the square of standart deviation coefficients a
|
|
- temp = c[0][0]*d[0][0] + c[1][1]*d[1][1] + c[2][2]*d[2][2]
|
|
|
|
|
|
+ temp = c[0][0] * d[0][0] + c[1][1] * d[1][1] + c[2][2] * d[2][2]
|
|
for i in range(3):
|
|
for i in range(3):
|
|
sigma2A.append(d[i][i] * S2 / temp)
|
|
sigma2A.append(d[i][i] * S2 / temp)
|
|
for i in range(len(data[4])): # recalculation of weight coefficients P
|
|
for i in range(len(data[4])): # recalculation of weight coefficients P
|
|
- data[4][i] = 1/(data[0][i]**2 * sigma2A[0] + sigma2A[1] + data[0][i]**2 * sigma2A[2] * (abs(data[3][i])**2))
|
|
|
|
|
|
+ data[4][i] = 1 / (
|
|
|
|
+ data[0][i] ** 2 * sigma2A[0] + sigma2A[1] + data[0][i] ** 2 * sigma2A[2] * (abs(data[3][i]) ** 2))
|
|
if error:
|
|
if error:
|
|
return abs(np.array(sigma2A))
|
|
return abs(np.array(sigma2A))
|
|
else:
|
|
else:
|
|
return data
|
|
return data
|
|
|
|
|
|
|
|
|
|
|
|
+def recalculating(data, a, c, d, n, printing=False):
|
|
|
|
+ for i in range(2, n):
|
|
|
|
+ data = recalculation_of_data(data, a, c, d)
|
|
|
|
+ a, c, d = solution(data)
|
|
|
|
+ Ql, diam, k, Q = q_factor(a)
|
|
|
|
+ sigma2A = recalculation_of_data(data, a, c, d, error=True)
|
|
|
|
+ sigmaQ0, sigmaQl = random_deviation(a, sigma2A, diam, k, Ql)
|
|
|
|
+ if printing:
|
|
|
|
+ print(f"Q = {Q} +- {sigmaQ0}, if i == {i}")
|
|
|
|
+
|
|
|
|
+
|
|
def random_deviation(a, sigma2A, diam, k, Ql):
|
|
def random_deviation(a, sigma2A, diam, k, Ql):
|
|
"""defines standart deviations of values"""
|
|
"""defines standart deviations of values"""
|
|
- sigmaQl = sigma2A[2]**0.5
|
|
|
|
- sigmaDiam = (sigma2A[0]/(abs(a[2])**2) + sigma2A[1] + abs(a[0]/a[2]/a[2])**2 * sigma2A[2])**0.5
|
|
|
|
- sigmaK = 2*sigmaDiam/((2-diam)**2)
|
|
|
|
- sigmaQ0 = ((1 + k)**2 * sigma2A[2] + Ql**2 * sigmaK**2)**0.5
|
|
|
|
|
|
+ sigmaQl = sigma2A[2] ** 0.5
|
|
|
|
+ sigmaDiam = (sigma2A[0] / (abs(a[2]) ** 2) + sigma2A[1] + abs(a[0] / a[2] / a[2]) ** 2 * sigma2A[2]) ** 0.5
|
|
|
|
+ sigmaK = 2 * sigmaDiam / ((2 - diam) ** 2)
|
|
|
|
+ sigmaQ0 = ((1 + k) ** 2 * sigma2A[2] + Ql ** 2 * sigmaK ** 2) ** 0.5
|
|
return sigmaQ0, sigmaQl
|
|
return sigmaQ0, sigmaQl
|
|
|
|
|
|
|
|
|
|
@@ -96,13 +110,7 @@ def apply(filename):
|
|
freq, re, im = open_file(filename)
|
|
freq, re, im = open_file(filename)
|
|
data = prepare_data(freq, re, im)
|
|
data = prepare_data(freq, re, im)
|
|
a, c, d = solution(data)
|
|
a, c, d = solution(data)
|
|
- for i in range(2, 10):
|
|
|
|
- data = recalculation_of_data(data, a, c, d)
|
|
|
|
- a, c, d = solution(data)
|
|
|
|
- Ql, diam, k, Q = q_factor(a)
|
|
|
|
- sigma2A = recalculation_of_data(data, a, c, d, error=True)
|
|
|
|
- sigmaQ0, sigmaQl = random_deviation(a, sigma2A, diam, k, Ql)
|
|
|
|
- print(f"Q = {Q} +- {sigmaQ0}, if i == {i}")
|
|
|
|
|
|
+ recalculating(data, a, c, d, 10, printing=True)
|
|
|
|
|
|
|
|
|
|
def fl_fitting(freq, re, im):
|
|
def fl_fitting(freq, re, im):
|
|
@@ -110,42 +118,35 @@ def fl_fitting(freq, re, im):
|
|
|
|
|
|
data, fl = prepare_data(freq, re, im)
|
|
data, fl = prepare_data(freq, re, im)
|
|
a, c, d = solution(data)
|
|
a, c, d = solution(data)
|
|
-
|
|
|
|
|
|
+ Ql, Q, sigmaQ0, sigmaQl = None, None, None, None
|
|
# Repeated curve fitting
|
|
# Repeated curve fitting
|
|
# 1.189 of Qfactor Matlab
|
|
# 1.189 of Qfactor Matlab
|
|
# fl2 = 0
|
|
# fl2 = 0
|
|
# g_d=0
|
|
# g_d=0
|
|
# g_c=0
|
|
# g_c=0
|
|
for x in range(0, 3):
|
|
for x in range(0, 3):
|
|
- g_c = (np.conj(a[2])*a[1]-a[0])/(np.conj(a[2])-a[2])
|
|
|
|
- g_d = a[0]/a[2]
|
|
|
|
- g_2 = 2*g_c-g_d
|
|
|
|
- dt = (a[1]-g_2)/(g_2*a[2]-a[0])
|
|
|
|
- fl2 = fl*(1 + np.real(dt)/2)
|
|
|
|
|
|
+ g_c = (np.conj(a[2]) * a[1] - a[0]) / (np.conj(a[2]) - a[2])
|
|
|
|
+ g_d = a[0] / a[2]
|
|
|
|
+ g_2 = 2 * g_c - g_d
|
|
|
|
+ dt = (a[1] - g_2) / (g_2 * a[2] - a[0])
|
|
|
|
+ fl2 = fl * (1 + np.real(dt) / 2)
|
|
data, fl = prepare_data(freq, re, im, fl2)
|
|
data, fl = prepare_data(freq, re, im, fl2)
|
|
a, c, d = solution(data)
|
|
a, c, d = solution(data)
|
|
-
|
|
|
|
- for i in range(2, 20):
|
|
|
|
- data = recalculation_of_data(data, a, c, d)
|
|
|
|
- a, c, d = solution(data)
|
|
|
|
-
|
|
|
|
- Ql, diam, k, Q = q_factor(a)
|
|
|
|
- sigma2A = recalculation_of_data(data, a, c, d, error=True)
|
|
|
|
- sigmaQ0, sigmaQl = random_deviation(a, sigma2A, diam, k, Ql)
|
|
|
|
|
|
+ recalculating(data, a, c, d, 20)
|
|
|
|
|
|
# taking into account coupling losses on page 69 of Qfactor Matlab
|
|
# taking into account coupling losses on page 69 of Qfactor Matlab
|
|
# to get results similar to example program
|
|
# to get results similar to example program
|
|
- if False:
|
|
|
|
- phi1=np.arctan(np.double(g_d.imag/g_d.real)) # 1.239
|
|
|
|
- phi2=np.arctan(np.double((g_c.imag-g_d.imag)/(g_c.real-g_d.real)))
|
|
|
|
- phi=-phi1+phi2
|
|
|
|
- d_s=(1-np.abs(g_d)**2)/(1-np.abs(g_d)*np.cos(phi))
|
|
|
|
- diam = abs(a[1] - a[0] / a[2])
|
|
|
|
- qk=1/(d_s/diam-1)
|
|
|
|
-
|
|
|
|
- sigma2A = recalculation_of_data(data, a, c, d, error=True)
|
|
|
|
- sigmaQ0 = random_deviation(a, sigma2A, diam, k, Ql)
|
|
|
|
- Q = Ql * (1 + qk) # Q-factor = result
|
|
|
|
- print(f"Q0 = {Q} +- {sigmaQ0}")
|
|
|
|
-
|
|
|
|
- return Q,sigmaQ0, Ql, sigmaQl,a
|
|
|
|
|
|
+ # if False:
|
|
|
|
+ # phi1=np.arctan(np.double(g_d.imag/g_d.real)) # 1.239
|
|
|
|
+ # phi2=np.arctan(np.double((g_c.imag-g_d.imag)/(g_c.real-g_d.real)))
|
|
|
|
+ # phi=-phi1+phi2
|
|
|
|
+ # d_s=(1-np.abs(g_d)**2)/(1-np.abs(g_d)*np.cos(phi))
|
|
|
|
+ # diam = abs(a[1] - a[0] / a[2])
|
|
|
|
+ # qk=1/(d_s/diam-1)
|
|
|
|
+ #
|
|
|
|
+ # sigma2A = recalculation_of_data(data, a, c, d, error=True)
|
|
|
|
+ # sigmaQ0 = random_deviation(a, sigma2A, diam, k, Ql)
|
|
|
|
+ # Q = Ql * (1 + qk) # Q-factor = result
|
|
|
|
+ # print(f"Q0 = {Q} +- {sigmaQ0}")
|
|
|
|
+
|
|
|
|
+ return Q, sigmaQ0, Ql, sigmaQl, a
|