フィルタの周波数特性¶
# モジュールの読み込み
import numpy as np
import pandas as pd
from scipy.special import comb
from scipy.special import eval_chebyt
import matplotlib.pyplot as plt
from matplotlib import style
import skrf as rf
import warnings
warnings.filterwarnings("ignore")
import scienceplots
plt.style.use(['science', 'notebook'])
#plt.rcParams['font.family'] = 'Times New Roman' # font familyの設定
plt.rcParams['font.family'] = 'serif' # font familyの設定
plt.rcParams['mathtext.fontset'] = 'cm' # math fontの設定
def dbdeg(z,iunwrap): # 電磁界の複素係数から絶対値，デジベル(dB)，偏角(deg)の計算
q = 180.0/np.pi
amp = np.abs(z)
db = 20.0*np.log10(amp)
rad = np.angle(z)
if iunwrap==1:
# rad = np.unwrap(rad, discont=np.pi)
rad = np.unwrap(rad)
deg = rad*q
return amp,db,deg
特性関数¶
動作伝送関数（Transducer Function） $H(s)$ は，特性関数（characteristic function）$K(s)$ を用いて（$s = j\omega$）， \begin{gather} |H(s)|^2= 1+|K(s)|^2 \end{gather} さらに，ポート1から見た反射係数を $\Gamma_1$ とすると， \begin{gather} |\Gamma_1|^2 = 1-\frac{1}{|H|^2} = \frac{|H|^2-1}{|H|^2} = \frac{|K|^2}{|H|^2} = \frac{|K|^2}{1-|K|^2} \end{gather}
最平坦特性¶
最平坦特性の特性関数 $K(s)$ は， \begin{gather} K(s) = s^N \end{gather} 遮断角周波数$\omega _c$で角周波数$\omega$を規格化して， 規格化角周波数$\Omega$を次式によって定義する． \begin{gather} s = j\Omega, \ \ \ \Omega = \frac{\omega}{\omega _c} \end{gather} これより，$|H|^2$は， \begin{gather} |H|^2 = 1 + (-s^2)^N = 1+ \Omega ^{2N} = 1 +\left( \frac{\omega}{\omega _c} \right) ^{2N} \end{gather} ただし，$\Omega =1$，$\omega = \omega _c$を遮断点と呼び， このとき，$|H(\Omega)|^2$は， $|H(1)|^2 = 1+1^{2N} = 2$である． また，$\Omega =0$のとき，$|H(0)|^2 = 1$． $N \to \infty $では，$0 \leq \Omega \leq 1$のとき $|H|^2= 1$で無損失伝送， $\Omega \geq 1$のとき$|H|^2=\infty$で損失が無限大（全く伝送しない）となり， 理想的な低域通過フィルタ（lowpass filter）の特性となる．
def maximally_flat_response(n,omega): # 最平坦特性の動作伝送関数 H, 反射係数の逆数
n2 = 2.0*n
h2 = 1.0+omega**n2
gamma2 = 1.0-1.0/h2
gamma2i = 1.0/gamma2
return h2, gamma2i
動作伝送関数の因数分解¶
実際の回路で実現できるように， $H(s)$を複素平面の左半平面内のゼロ点の値を用いて因数分解すると， \begin{gather} H(s) = (s-j\epsilon ) (s-j\epsilon ^3) \cdots (s-j\epsilon ^{2N-1}) = \prod _{n=1}^N \left( s-j\epsilon ^{2n-1} \right) = \prod _{n=1}^N \left( s-s_n \right) \end{gather} ここで， \begin{gather} s_n = j \epsilon ^{2n-1} = j e^{j\frac{2n-1}{2N} \pi} = e^{j\frac{2n-1+N}{2N} \pi} \end{gather} このとき，反射係数$\Gamma(s)$は， \begin{gather} \Gamma(s) = \frac{\pm K(s)}{H(s)} = \pm \frac{s^N}{\displaystyle{\prod _{n=1}^N \left( s-s_n \right)}} \end{gather} ただし，上式の符号は，この時点では定まらないので，プログラムではプラス符号を計算する．後述する等リプル特性の場合も同様である．
def maximally_flat_response_h(n,omega): # 因数分解した最平坦特性の動作伝送関数 H, 反射係数，sn
s = 1.0j*omega
h = 1.0
ssn = np.empty(0)
for i in range(1,n+1):
th = (2.0*i-1)/2.0/n*np.pi
sn = 1.0j*np.exp(1.0j*th)
ssn = np.append(ssn, sn)
h = h*(s-sn)
gammai = h/s**n
return h, gammai, ssn
def maximally_flat_ssn(n,omega):
ssn = np.empty(0)
for i in range(1,n+1):
th = (2.0*i-1)/2.0/n*np.pi
sn = 1.0j*np.exp(1.0j*th)
ssn = np.append(ssn, sn)
return ssn
最平坦特性を持つ低域通過梯子型回路の規格化素子値¶
低域通過梯子型回路は双対的な次のような回路で表され，$N$次の最平坦特性を有する梯子型リアクタンス回路を合成すると，規格化素子値$g_k$は，次のようになる． \begin{align} g_0 &= 1\\ g_k &= 2 \sin \left( \frac{(2k-1) \pi }{2N} \right) \ \ \ (k=1,2, \cdots ,N)\\ g_{N+1} &= 1 \end{align} ただし，$g_0$は入力の負荷抵抗，$g_{N+1}$は終端負荷の抵抗を示す．
def element_values_maximally_flat(n): # 最平坦特性を持つ低域通過梯子型回路の規格化素子値
gg = np.empty(0)
gg = np.append(gg, 1.0)
for k in range(1,n+1):
th = (2.0*k-1.0)*np.pi/2.0/n
g = 2.0*np.sin(th)
gg = np.append(gg, g)
gg = np.append(gg, 1.0)
return gg
等リプル特性¶
特性関数$K(s)$を，$s = j \Omega$ より， \begin{gather} K(s) = K(j \Omega) \equiv \varepsilon T_N(\Omega) \end{gather} ここで，$T_N(\Omega)$は$N$次のチェビシェフ多項式（Chebyshev polynomial）を示し， \begin{gather} T_N (\Omega) = \left\{ \begin {array}{ll} (-1)^N \cosh \big( N \cosh ^{-1} |\Omega| \big) & (\Omega < -1) \\ \cos \big( N \cos^{-1} \Omega \big) & (|\Omega| \le 1) \\ \cosh \big( N \cosh ^{-1} \Omega \big) & (\Omega > 1) \end{array} \right. \end{gather} また，$\varepsilon$ は通過域のリプルの大きさを決めるパラメータ（実数）である．これより， \begin{align} |K(s)|^2 &= K(s) K(s)^*=K(s) K(-s) =K(j\Omega) K(-j\Omega) \nonumber \\ &= \varepsilon T_N(\Omega) \cdot \varepsilon T_N(-\Omega) = \varepsilon^2 T_N^2(\Omega) \\ |H(s)|^2 &= 1+|K(s)|^2 = 1+\varepsilon^2 T_N^2 (\Omega) \end{align} ここで，遮断点$\Omega =1$のとき，$|H|^2$は，$T_N(1)=1$より， \begin{gather} |H|^2 \Big| _{\Omega=1} = 1+\varepsilon^2 T_N^2 (1) = 1+\varepsilon^2 \end{gather} 通過域のリプルの最大値は， \begin{gather} |H|_{max} = \sqrt{1+\varepsilon^2} \end{gather} デジベル値は，$L_{Ar}=10 \log_{10} (1+\varepsilon^2)$ [dB]． 逆に，$\varepsilon$は， \begin{gather} \varepsilon = \sqrt{10^{\frac{L_{Ar}}{10}}-1} \end{gather}
def chebyshev_response(n,rwdb,omega):
k2 = 10.0**(rwdb/10.0)-1.0
che = eval_chebyt(n, omega)
h2 = 1.0+k2*che**2
gamma2 = 1.0-1.0/h2
gamma2i = 1.0/gamma2
return h2, gamma2i
動作伝送関数の因数分解¶
動作伝送関数$H(s)$，さらには反射係数$\Gamma(s)$は， \begin{gather} H(s) = \pm \varepsilon 2^{N-1} (-j)^N \prod _{n=1}^N \left( s-s_n \right) \\ \Gamma(s) %= \frac{K(s)}{H(s)} = \pm \frac{T_N (\Omega)}{\displaystyle{2^{N-1} (-j)^N \prod _{n=1}^N \left( s-s_n \right)}} \end{gather} ここで，$s_n$は， \begin{align} s_n &=\sigma_n + j \omega_n \\ \sigma_n &= - \sin \frac{(2n+1)\pi}{2N} \sinh a \\ \omega_n &= \cos \frac{(2n+1)\pi}{2N} \cosh a \\ a &= \frac{1}{N} \sinh ^{-1} \frac{1}{\varepsilon} \end{align} このとき，次の関係が成り立つ． \begin{gather} \sin^2 \frac{(2n+1)\pi}{2N} + \cos^2 \frac{(2n+1)\pi}{2N} = \frac{\sigma_n^2}{\sinh ^2a} + \frac{\omega_n^2}{\cosh ^2a} = 1 \end{gather} これより，$s_n$は複素平面のだ円（直交するだ円の軸は実軸と虚軸）上にあることがわかる．
def chebyshev_response_h_(n,rwdb,omega): # 因数分解した等リプル特性の動作伝送関数 H, 反射係数，sn
s = 1.0j*omega
k2 = 10.0**(rwdb/10.0)-1.0
k = np.sqrt(k2)
ki = 1.0/k
a = np.sqrt(1.0+ki**2)+ki
an = a**(1.0/n)
ani = 1.0/an
sh = (an-ani)/2.0
ch = (an+ani)/2.0
# print("sinh a =",sh, ", cosh a=",ch)
h = 2.0**(n-1)*k
ssn = np.empty(0)
for i in range(1,n+1):
th = (2.0*i-1.0)/2.0/n*np.pi
re = -sh*np.sin(th)
im = ch*np.cos(th)
sn = re+1.0j*im
ssn = np.append(ssn, sn)
h = h*(s-sn)
che = eval_chebyt(n, omega)
fn = 1.0j**n*che
gammai = h/k/fn
return h, gammai, ssn, sh, ch
def chebyshev_ssn(n,rwdb,omega):
k2 = 10.0**(rwdb/10.0)-1.0
k = np.sqrt(k2)
ki = 1.0/k
a = np.sqrt(1.0+ki**2)+ki
an = a**(1.0/n)
ani = 1.0/an
sh = (an-ani)/2.0
ch = (an+ani)/2.0
ssn = np.empty(0)
for i in range(1,n+1):
th = (2.0*i-1.0)/2.0/n*np.pi
re = -sh*np.sin(th)
im = ch*np.cos(th)
sn = re+1.0j*im
ssn = np.append(ssn, sn)
return ssn, sh, ch
等リプル特性を持つ低域通過梯子型回路の規格化素子値¶
$N$次のチェビシェフ特性を有する梯子型リアクタンス回路を合成すると，規格化素子値$g_k$は， \begin{align} g_1 &= \frac{2 u_1}{\eta} \\ g_k &= \frac{4u_{k-1} u_k}{v_{k-1} g_{k-1}} \ \ \ (k=2,3, \cdots ,N) \\ g_{N+1} &= \left\{ \begin {array}{cc} \displaystyle{1} & (N:\mbox{奇数}) \\ \displaystyle{\coth ^2 \frac{\zeta}{4}} & (N:\mbox{偶数}) \end{array} \right. \end{align} ただし，$g_0 =1$ である．また，$g_{N+1}$は終端負荷の抵抗（最終段が直列素子の場合）あるいはコンダクタンス（最終段が並列素子の場合）を示す．ここで， \begin{align} \zeta &= \ln \left[ \coth \left( \frac{L_{Ar}}{17.37} \right) \right] \\ \eta &= \sinh \left( \frac{\zeta}{2N} \right) \\ u_k &= \sin \left( \frac{(2k-1)\pi}{2N} \right) \ \ \ (k=1,2, \cdots ,N) \\ v_k &= \eta ^2 + \sin ^2 \left( \frac{k\pi}{N} \right) \ \ \ (k=1,2, \cdots ,N) \end{align} ただし，$L_{Ar}$はリプルのピーク値[dB]を示す．
def element_values_equal_ripple(nn,rwdb):
ch = 1.0/np.tanh(rwdb/17.37)
beta = np.log(ch)
gamma = np.sinh(beta/2.0/nn)
gg = np.empty(0)
gg = np.append(gg, 1.0)
for i in range(1,nn+1):
th = (2.0*i-1.0)*np.pi/2.0/nn
ai = np.sin(th)
bi = gamma**2+(np.sin(i*np.pi/nn))**2
if i == 1:
gi = 2.0*ai/gamma
else:
gi = 4.0*aim*ai/bim/gim
aim = ai
bim = bi
gim = gi
gg = np.append(gg, gi)
if nn % 2 == 1:
gn1 = 1.0
else:
gn1 = (1.0/np.tanh(beta/4.0))**2
gg = np.append(gg, gn1)
return gg
入力データ¶
ic = 1 # maximally flat
# ic, rwdb = 2, 0.04321 # equal_ripple, リプル [dB]
# ic, rwdb = 2, -10.0*np.log10(1.0-10.0**(-20.0/10.0)) # equal_ripple, リプル [dB]
if ic == 2:
print("RW=",rwdb,"[dB], s11=",-10.0*np.log10(1.0-10.0**(-rwdb/10.0)),"[dB]")
ip = 1 # low pass
# ip = 2 # high pass
# ip, rbw = 3, 0.3 # band pass
# ip, rbw = 4, 0.3 # band stop
# 段数
nn = np.array([3, 5, 7])
# nn = np.arange(2, 10, 1)
# nn = np.arange(3, 8, 2)
# nn = np.array([2,3,4,5])
n_n = len(nn)
print("nn=",nn)
nn= [3 5 7]
n_ome = 501
omega = np.linspace(0.01, 0.1, n_ome)
omega = np.append(omega, np.linspace(0.1, 1.0, n_ome))
omega = np.append(omega, np.linspace(1.0, 3.0, n_ome))
omega = np.append(omega,np.linspace(3.0, 10.0, n_ome))
n_w = len(omega)
print("n_w=",n_w)
n_w= 2004
周波数変換¶
if ip == 1: # 低域通過
omega_t = 1.0*omega
if ip == 2: # 低域通過から高域通過への周波数変換
omega_t = -1.0/omega
elif ip == 3: # 低域通過から帯域通過への周波数変換
omega_t = (omega-1.0/omega)/rbw
elif ip == 4: # 低域通過から帯域阻止への周波数変換
omega_t = (omega-1.0/omega)/rbw
omega_t = np.where(omega_t==0.0, 1.0e-10, omega_t)
omega_t = -1.0/omega_t
hh2 = np.empty((0,n_w))
gg2 = np.empty((0,n_w))
hh = np.empty((0,n_w))
ggi = np.empty((0,n_w))
for n in nn:
if ic == 1:
gk = element_values_maximally_flat(n) # 特性関数K(s)=s**n
h2, g2 = maximally_flat_response(n,omega_t)
h, gammai, ssn = maximally_flat_response_h(n,omega_t)
elif ic == 2:
gk = element_values_equal_ripple(n,rwdb)
h2, g2 = chebyshev_response(n,rwdb,omega_t)
h, gammai, ssn, sh1, ch1 = chebyshev_response_h_(n,rwdb,omega_t)
n_ssn = len(ssn)
print(f'{n:2d}',end=' ')
[print(f"{gk[j]:6.4f}",end=' ') for j in range(n+2)]
print('')
hh2 = np.append(hh2, [h2], axis=0)
gg2 = np.append(gg2, [g2], axis=0)
hh = np.append(hh, [h], axis=0)
ggi = np.append(ggi, [gammai], axis=0)
hh2_amp = np.abs(hh2)
hh2_db = 10.0*np.log10(hh2_amp)
gg2_amp = np.abs(gg2)
gg2_db = 10.0*np.log10(gg2_amp)
hh_amp,hh_db,hh_deg = dbdeg(hh,0)
ggi_amp,ggi_db,ggi_deg = dbdeg(ggi,0)
gg = 1.0/ggi
3 1.0000 1.0000 2.0000 1.0000 1.0000 5 1.0000 0.6180 1.6180 2.0000 1.6180 0.6180 1.0000 7 1.0000 0.4450 1.2470 1.8019 2.0000 1.8019 1.2470 0.4450 1.0000
挿入損失 [dB]¶
df = pd.DataFrame({
f'$N={N}$': hh_db[i, ::100] for i, N in enumerate(nn)
})
df.index = omega[::100]
df.index.name = '$\Omega$'
df.columns.name = '$|H|$ [dB]'
df.style.background_gradient(cmap='rainbow')
| $|H|$ [dB] | $N=3$ | $N=5$ | $N=7$ |
|---|---|---|---|
| $\Omega$ | |||
| 0.010000 | 0.000000 | -0.000000 | -0.000000 |
| 0.028000 | 0.000000 | 0.000000 | -0.000000 |
| 0.046000 | 0.000000 | 0.000000 | -0.000000 |
| 0.064000 | 0.000000 | 0.000000 | 0.000000 |
| 0.082000 | 0.000001 | 0.000000 | 0.000000 |
| 0.100000 | 0.000004 | 0.000000 | 0.000000 |
| 0.278200 | 0.002013 | 0.000012 | 0.000000 |
| 0.458200 | 0.040005 | 0.001771 | 0.000078 |
| 0.638200 | 0.283955 | 0.048409 | 0.008068 |
| 0.818200 | 1.139516 | 0.547894 | 0.254126 |
| 0.998200 | 2.986890 | 2.971354 | 2.955875 |
| 1.392000 | 9.177701 | 14.520080 | 20.151636 |
| 1.792000 | 15.329486 | 25.346500 | 35.468554 |
| 2.192000 | 20.489408 | 34.085750 | 47.717750 |
| 2.592000 | 24.832397 | 41.363817 | 57.908907 |
| 2.992000 | 28.563745 | 47.596234 | 66.634623 |
| 4.358000 | 38.357867 | 63.928724 | 89.500212 |
| 5.758000 | 45.616419 | 76.027166 | 106.438032 |
| 7.158000 | 51.287534 | 85.479169 | 119.670837 |
| 8.558000 | 55.942348 | 93.237228 | 130.532120 |
| 9.958000 | 59.890332 | 99.817212 | 139.744097 |
透過係数の位相特性（主値）[deg]¶
df = pd.DataFrame({
f'$N={N}$': -hh_deg[i, ::100] for i, N in enumerate(nn)
})
df.index = omega[::100]
df.index.name = '$\Omega$'
df.columns.name = 'Arg $(1/H)$ [deg]'
df.style.background_gradient(cmap='coolwarm')
| Arg $(1/H)$ [deg] | $N=3$ | $N=5$ | $N=7$ |
|---|---|---|---|
| $\Omega$ | |||
| 0.010000 | -1.145935 | -1.854154 | -2.574880 |
| 0.028000 | -3.208983 | -5.192083 | -7.210250 |
| 0.046000 | -5.273075 | -8.531300 | -11.847289 |
| 0.064000 | -7.338891 | -11.872635 | -16.487077 |
| 0.082000 | -9.407123 | -15.216928 | -21.130698 |
| 0.100000 | -11.478482 | -18.565027 | -25.779249 |
| 0.278200 | -32.326600 | -52.110794 | -72.314170 |
| 0.458200 | -54.729324 | -87.516479 | -121.226344 |
| 0.638200 | -79.662947 | -126.334273 | -174.206298 |
| 0.818200 | -107.291566 | -172.099149 | 123.530148 |
| 0.998200 | -134.741937 | 135.513247 | 45.795625 |
| 1.392000 | -178.271451 | 55.626772 | -69.729422 |
| 1.792000 | 158.184256 | 18.436889 | -120.107272 |
| 2.192000 | 144.469126 | -2.889910 | -149.332770 |
| 2.592000 | 135.480031 | -16.997354 | -168.781919 |
| 2.992000 | 129.099932 | -27.095556 | 177.258705 |
| 4.358000 | 116.539297 | -47.161673 | 149.461742 |
| 5.758000 | 110.004863 | -57.673534 | 134.880272 |
| 7.158000 | 106.062157 | -64.032100 | 126.055832 |
| 8.558000 | 103.420958 | -68.296620 | 120.136192 |
| 9.958000 | 101.527061 | -71.356469 | 115.888241 |
反射損失 [dB]¶
df = pd.DataFrame({
f'$N={N}$': ggi_db[i, ::100] for i, N in enumerate(nn)
})
df.index = omega[::100]
df.index.name = '$\Omega$'
df.columns.name = '$1/|\Gamma_1|$ [dB]'
df.style.background_gradient(cmap='rainbow')
| $1/|\Gamma_1|$ [dB] | $N=3$ | $N=5$ | $N=7$ |
|---|---|---|---|
| $\Omega$ | |||
| 0.010000 | 120.000000 | 200.000000 | 280.000000 |
| 0.028000 | 93.170518 | 155.284197 | 217.397876 |
| 0.046000 | 80.234530 | 133.724217 | 187.213904 |
| 0.064000 | 71.629202 | 119.382003 | 167.134804 |
| 0.082000 | 65.171170 | 108.618615 | 152.066061 |
| 0.100000 | 60.000004 | 100.000000 | 140.000000 |
| 0.278200 | 33.340585 | 55.564299 | 77.790002 |
| 0.458200 | 20.376700 | 33.896263 | 47.452366 |
| 0.638200 | 11.986547 | 19.552729 | 27.314116 |
| 0.818200 | 6.367947 | 9.261946 | 12.453799 |
| 0.998200 | 3.033836 | 3.049598 | 3.065416 |
| 1.392000 | 0.559347 | 0.156157 | 0.042143 |
| 1.792000 | 0.129205 | 0.012699 | 0.001233 |
| 2.192000 | 0.038975 | 0.001695 | 0.000073 |
| 2.592000 | 0.014297 | 0.000317 | 0.000007 |
| 2.992000 | 0.006049 | 0.000076 | 0.000001 |
| 4.358000 | 0.000634 | 0.000002 | 0.000000 |
| 5.758000 | 0.000119 | 0.000000 | 0.000000 |
| 7.158000 | 0.000032 | 0.000000 | 0.000000 |
| 8.558000 | 0.000011 | 0.000000 | 0.000000 |
| 9.958000 | 0.000004 | 0.000000 | 0.000000 |
fig1 = plt.figure(figsize=(6, 10)) # グラフ領域の作成
plt.subplot(3, 1, 1) # 上側の図
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.01, 10.0) # x軸範囲の設定
plt.xticks(fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylim(0.0, 50.0) # y軸範囲の設定
plt.yticks(np.arange(0.0, 60.0, step=10.0), fontsize=11)
plt.ylabel("Insertion loss [dB]", fontsize=11) # y軸のラベル設定
plt.xscale('log') # log スケール
for i in range(n_n):
# plt.plot(omega, hh2_db[i,:], label=f"$N=${nn[i]:2d}")
# plt.plot(omega, gg2_db[i,:], label=f"$1/|\Gamma|^2: N=${nn[i]:2d}")
plt.plot(omega, hh_db[i,:], label=f"$|H|: N=${nn[i]:2d}")
# plt.plot(omega, ggi_db[i,:], label=f"$1/|\Gamma|^2: N=${nn[i]:2d}")
plt.legend(ncol=1, loc='upper left', fancybox=False, frameon = True, fontsize=9)
plt.subplot(3, 1, 2) # 下側の図
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.01, 10.0) # x軸範囲の設定
plt.xticks(fontsize=11)
plt.ylim(-180.0, 180.0) # y軸範囲の設定
plt.yticks(np.arange(-180.0, 270, step=90.0), fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylabel("Arg $(1/H)$ [deg]", fontsize=11) # y軸のラベル設定
plt.xscale('log') # log スケール
for i in range(n_n):
plt.plot(omega, -hh_deg[i,:], label=f"$N=${nn[i]:2d}", marker = '.', linewidth=0)
plt.legend(ncol=1, loc='lower left', fancybox=False, frameon = True, fontsize=9)
plt.subplot(3, 1, 3)
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.01, 10.0) # x軸範囲の設定
plt.xticks(fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylim(0.0, 50.0) # y軸範囲の設定
plt.yticks(np.arange(0.0, 60.0, step=10.0), fontsize=11)
plt.ylabel("Reflection loss [dB]", fontsize=11) # y軸のラベル設定
plt.xscale('log') # log スケール
for i in range(n_n):
# plt.plot(omega, hh2_db[i,:], label=f"$N=${nn[i]:2d}")
# plt.plot(omega, gg2_db[i,:], label=f"$1/|\Gamma|^2: N=${nn[i]:2d}")
# plt.plot(omega, hh_db[i,:], label=f"$|H|^2: N=${nn[i]:2d}")
plt.plot(omega, ggi_db[i,:], label=f"$1/|\Gamma|^2: N=${nn[i]:2d}")
plt.legend(ncol=1, loc='upper left', fancybox=False, frameon = True, fontsize=9)
plt.tight_layout()
plt.show()
挿入損失，位相特性，反射損失¶
hh_amp,hh_db,hh_deg = dbdeg(hh,1)
fig2 = plt.figure(figsize=(9, 6)) # グラフ領域の作成
plt.subplot(2, 2, 1) # 左上側の図
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.0, 2.0) # x軸範囲の設定
plt.xticks(np.arange(0.0, 2.5, step=0.5), fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylim(0.0, 3.0) # y軸範囲の設定
plt.yticks(np.arange(0.0, 3.5, step=0.5), fontsize=11)
plt.ylabel("$|H|^2$ [dB]", fontsize=11) # y軸のラベル設定
for i in range(n_n):
plt.plot(omega, hh2_db[i,:], label=f"$N=${nn[i]:2d}")
# plt.plot(omega, gg2_amp[i,:], label=f"$1/|\Gamma|^2: N=${nn[i]:2d}")
plt.legend(ncol=1, loc='upper right', fancybox=False, frameon = True, fontsize=11)
plt.subplot(2, 2, 2) # 右上の図
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.0, 2.0) # x軸範囲の設定
plt.xticks(np.arange(0.0, 2.5, step=0.5), fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylim(0.0, 40.0) # y軸範囲の設定
plt.yticks(np.arange(0.0, 50.0, step=10.0), fontsize=11)
plt.ylabel("$|H|^2$ [dB]", fontsize=8) # y軸のラベル設定
for i in range(n_n):
plt.plot(omega, hh2_db[i,:], label=f"$N=${nn[i]:2d}")
# plt.plot(omega, gg2_db[i,:], label=f"$1/|\Gamma|^2: N=${nn[i]:2d}")
plt.legend(ncol=1, loc='upper left', fancybox=False, frameon = True, fontsize=11)
plt.tight_layout()
plt.subplot(2, 2, 3) # 左下側の図
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.0, 2.0) # x軸範囲の設定
plt.xticks(np.arange(0.0, 2.5, step=0.5), fontsize=11)
plt.ylim(-1440.0, 180.0) # y軸範囲の設定
plt.yticks(np.arange(-1440.0, 360, step=180.0), fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylabel("Arg $(1/H)$ [deg]", fontsize=11) # y軸のラベル設定
for i in range(n_n):
plt.plot(omega, -hh_deg[i,:], label=f"$N=${nn[i]:2d}", marker = '.', linewidth=0)
plt.legend(ncol=1, loc='lower left', fancybox=False, frameon = True, fontsize=11)
plt.tight_layout()
plt.subplot(2, 2, 4) # 右下の図
plt.grid(color = "gray", linestyle="--")
plt.minorticks_on() # #補助目盛りをつける
plt.xlim(0.0, 2.0) # x軸範囲の設定
plt.xticks(np.arange(0.0, 2.5, step=0.5), fontsize=11)
plt.xlabel("$\Omega$", fontsize=11) # x軸のラベル設定
plt.ylim(0.0, 40.0) # y軸範囲の設定
plt.yticks(np.arange(0.0, 50.0, step=10.0), fontsize=11)
plt.ylabel("$1/|\Gamma|^2$ [dB]", fontsize=11) # y軸のラベル設定
for i in range(n_n):
plt.plot(omega, gg2_db[i,:], label=f"$N=${nn[i]:2d}")
plt.legend(ncol=1, loc='upper right', fancybox=False, frameon = True, fontsize=11)
plt.tight_layout()
plt.show()
$s_n$の位置¶
n = nn[1]
n
5
if ic == 1:
ssn = maximally_flat_ssn(n,omega_t)
elif ic == 2:
ssn, sh1, ch1 = chebyshev_ssn(n,rwdb,omega_t)
if ic==1:
n_th = 361
th = np.linspace(0.0, 2.0*np.pi, n_th)
cir = np.exp(1.0j*th)
s_param = cir.reshape(-1, 1, 1) # (n_freqs, 1, 1)
p5 = rf.Network(f=np.ones(n_th), s=s_param, z0=np.ones(n_th))
elif ic==2:
n_xp = 201
x_daen = np.linspace(-sh1, sh1, n_xp)
y_daen = np.sqrt(1.0-x_daen**2/sh1**2)*ch1
xy_daen1 = x_daen + 1.0j*y_daen
x_daen = np.linspace(sh1, -sh1, n_xp)
y_daen = np.sqrt(1.0-x_daen**2/sh1**2)*ch1
xy_daen2 = x_daen - 1.0j*y_daen
daen = np.concatenate((xy_daen1, xy_daen2)) # 上下の半楕円を連結
s_param = daen.reshape(-1, 1, 1) # (n_freqs, 1, 1)
p5 = rf.Network(f=np.ones(len(daen)), s=s_param, z0=np.ones(len(daen)))
with plt.style.context('seaborn-v0_8-ticks'):
fig3, ax = plt.subplots(figsize=(6, 6))
# 軌跡（p5）を描画
s_data = p5.s[:, 0, 0] # shape: (n_freqs,)
ax.plot(s_data.real, s_data.imag, linewidth=1, label='') # 凡例付き
# 離散点（ssn[i]）をプロットして凡例を設定
for i in range(len(ssn)):
point = ssn[i]
i1 = i+1
ax.plot(point.real, point.imag,
marker='o', markersize=7, linewidth=0,
label=f'$s_{i1}$') # 任意のラベル
ax.set_aspect('equal') # 円が歪まないように
ax.set_xlabel('Real Part', fontsize=11)
ax.set_ylabel('Imaginary Part', fontsize=11)
ax.tick_params(labelsize=11)
ax.set_xlim([-1.6, 1.6])
ax.set_ylim([-1.6, 1.6])
ax.grid(color="gray", linestyle="--")
ax.legend(ncol=1, loc='best', fancybox=False, frameon = True, fontsize=13) # 凡例を表示
plt.show()
スミス図表（反射係数）¶
n_w2 = n_w
fig4 = plt.figure(figsize=(9, 9))
ax = plt.subplot(1, 1, 1)
rf.plotting.smith(ax=ax, chart_type='zy', draw_labels=True) # smith座標系の準備
colors = plt.cm.tab10.colors # 10色までの色分け
labels = [] # 凡例ラベル格納用
for i in range(n_n):
s_param = gg[i, :].reshape(-1, 1, 1)
p3 = rf.Network(f=omega, s=s_param[:n_w2], z0=np.ones(n_w2))
data = p3.s[:, 0, 0] # 1ポートの反射係数（complex）
line, = ax.plot(data.real, data.imag, linewidth=2, color=colors[i % 10]) # プロット
line.set_label(f'$N= {nn[i]}$') # ラベル設定
plt.legend(ncol=1, loc='best', fancybox=False, frameon = True, fontsize=12) # 凡例を表示
<matplotlib.legend.Legend at 0x31bacad10>
if ic == 1:
fig1.savefig('maximally_flat_amp_ph_n_.pdf')
fig2.savefig('maximally_flat_loss_n_.pdf')
fig3.savefig('maximally_flat_sn.pdf')
fig4.savefig('maximally_flat_smith_chart_n_.pdf')
elif ic ==2:
fig1.savefig('equal_ripple_amp_ph_n_.pdf')
fig2.savefig('equal_ripple_loss_n_.pdf')
fig3.savefig('equal_ripple_sn.pdf')
fig4.savefig('equal_ripple_smith_chart_n_.pdf')