同調法はopen/short 法の高域側限界を決める 伝送線路の共振を積極的に利用します。 共振周波数間隔と試料長から次式
v = 2 * l * Δf (1) ここに、 v = 伝送線路中の電磁波の位相速度 (m/s) l = 伝送線路の長さ (m) Δf = 伝送線路の入力インピーダンスが極小になるような共振周波数間隔 (Hz)で電磁波の速度がわかりますから、
Vr = v / c (2) ここに、 Vr = 速度係数 (0 < Vr <= 1) c = 真空中の光速 (2.99792458e8 m/s .. 定義値)で速度係数が得られ、必要なら
β = ω / v + 2 * n * π (3) ここに、 β = 位相定数 (rad/m) ω = 角周波数 (rad/s) = 2 * π * f f = 周波数 (Hz) n = 整数 (0, 1, 2, ..)で位相定数も決められますが、長さに依存しない速度係数のほうが使いやすいです。
一般に周波数の測定精度は極めて高いですから、 速度係数の測定誤差は長さの測定精度で決まります。
次ぎに、特性インピーダンスは高周波に於ける精度の高い近似式、
Z0 = 1 / (v * C) (4) ここに、 Z0 = 特性インピーダンス .. 純抵抗に極めて近い (Ω) C = 伝送線路のキャパシタンス (F/m)を使って、キャパシタンスがわかれば Z0 が得られます。 高周波の C 測定は無理ですが、 高周波用伝送線路に使われる誘電材料なら低周波の値と変わらないので、 高周波用インピーダンス測定機で使える最低周波数は 1 MHz 程度が普通ですから、 長さ 1 m の試料ならポリエチレン絶縁の同軸ケーブルでも、 最低の共振周波数 50 MHz 程度近くまでは伝送線路補正を使った C 測定が可能です。
最後に、減衰定数については、入力インピーダンスが極小になる共振周波数で、
α = atanh(Zmin / Z0) / l ここに、 α = 減衰定数 (neper/m) Zmin = 入力インピーダンスの極小値 (Ω) .. 純抵抗 l = 試料の長さ (m)で決められます。
かくて、高周波用インピーダンス測定機を使えば、 終端開放試料を一度接続するだけで、 試料の二次定数を測定することができることになります。 (注１)
以下、前期の方針に基づく測定のシミュレーションを行うプログラム例です。 長さ 1 m の終端開放 3C-2V 程度のポリエチレン絶縁同軸ケーブルの 入力インピーダンスを 1 MHz から 1 GHz まで測定できる インピーダンス測定機で測定する状況を想定しています。
キャパシタンスの決定は１変数の回帰分析を使っていますが、 回帰直線の傾きが 0 になるとき分散が最小になりますから、 短に分散を最小にする平均値を求めれば、 それがキャパシタンスになります。
誤差要因の追加法は既に説明しましたので、 実験してみると、 open/short 法と比べて、 かなり精度が高いことが実感できます。 共振法シミュレーションでも実験できるようにしておきました。
/usr/bin/awk '
BEGIN {
_Z0 = 75.0
_a = 1.373e-6
_b = 8.385e-12
_v = 1 / sqrt(2.3) # ポリエチレン絶縁
_R = 0.0914
_C = 69e-12
l = 1 # 試料長 (m)
PI = atan2(0, -1)
Vr = _v
# 1/4 波長共振周波数探索
f0o = find_f0o(l, Vr) # 4.94193e+07
print "f0o =", f0o >"/dev/stderr" # 1/4 波長共振周波数
# (2*n+1)/4 (n = 1, 2, ..) 共振周波数探索
K = 0
Fr[K] = f0o # 1/4 波長周波数
Alpha[K++] = x # 極小値 (純抵抗) .. 減衰定数計算の準備
for (i = 1; i < 20; i++) {
f = f0o * (2 * i + 1)
if (1e9 < f)
break
f = find_next(f - f0o * 0.20, f + f0o * 0.20)
print "#", f, f0, f - f0; f0 = f
Fr[K] = f
Alpha[K++] = x
}
df0 = 9.91764e+07 # 真値(参考)
df = (Fr[K - 1] - Fr[0]) / (K - 1)
print "df =", df, df / df0 >"/dev/stderr"
Vr = 2 * df * l / 2.99792458e8 # 速度係数
print "Vr =", Vr >"/dev/stderr"
# キャパシタンス測定
M = int(int(Fr[0] / 1e6 * 0.80)) # 最初の共振点直前まで
N = 0
for (i = 1; i <= M; i++) { # 1 MHz 刻みで測定
f = i * 1e6
fm[N] = f
meas_c(f) # Zin 測定
rm[N] = x; xm[N] = y
N++
}
cap(N) # 伝送線路補正
print "C =", C >"/dev/stderr"
# 特性インピーダンスは高周波の速度係数を使う
Z0 = 1 / C / Vr / 2.99792458e8
print "Z0 =", Z0 >"/dev/stderr"
print "Vr =", Vr, "C =", C, "Z0 =", Z0 >"/dev/stderr"
# 減衰定数を決定
print "x Freq. (Hz)"
print "y Att. (neper/m)"
printf("s Z0 = %.1f Ohm, Vr = %.3f\n", Z0, Vr)
# Alpha
for (i = 0; i < K; i++) {
Alpha[i] = atanh(Alpha[i] / Z0)
print Fr[i], Alpha[i]
}
# 減衰定数回帰式 Aloha(f) = b1 * sqrt(f) の b1 を求める
for (i = 0; i < K; i++) {
Y[i] = Alpha[i]
X[i] = sqrt(Fr[i])
}
reg0(K) # 回帰直線の係数を求める .. reg1(K) でもよい
print "b0 =", b0, "b1 =", b1 >"/dev/stderr"
print "Alpfa(f) =", b1, "*sqrt(f)" >"/dev/stderr"
# exit 0
# 回帰式 Aloha(f) = b0 + b1 * sqrt(f) + b2 がほしい場合
for (i = 0; i < K; i++) {
Y[i] = Alpha[i]
X1[i] = sqrt(Fr[i])
X2[i] = Fr[i]
}
reg2(K)
print "b0 =", b0, "b1 =", b1, "b2 =", b2 >"/dev/stderr"
print "Alpfa(f) =", b0, "+", b1, "* sqrt(f) +", b2, "* f" >"/dev/stderr"
exit 0
}
# キャパシタンス決定
function cap(M, i, f1, bl, cm, a, b, c, d, p, y1, y2) {
# 速度係数の決定
# (β*l)/tan(β*l) 補正したキャパシタンスが一定になる β を求める
# 黄金分割法で回帰直線の分散の極小値を見付ける
a = 0.45 # 0.45 <= vr <= 1 の範囲を探索
b = 1.00
p = (sqrt(5) - 1) * 0.5
c = b - p * (b - a)
d = a + p * (b - a)
y1 = cap2(c, M)
y2 = cap2(d, M)
for (i = 0; i < 100; i++) {
if ((d - c) / (c + d) < 1e-3)
break
if (y1 < y2) {
b = d
d = c
y2 = y1
c = b - p * (b - a)
y1 = cap2(c, M)
}
else if (y1 > y2) {
a = c
c = d
y1 = y2
d = a + p * (b - a)
y2 = cap2(d, M)
}
}
vr = (c + d) * 0.5
# 低域の速度係数は高域よりわずかに小さい
# 単位長あたりキャパシタンス
a = b = 0
for (i = 0; i < M; i++) {
w = 2 * PI * fm[i]
cm = -xm[i] / w / (rm[i] * rm[i] + xm[i] * xm[i])
bl = w / vr / 2.99792458e8 * l
c = cm / sin(bl) * cos(bl) * bl
a += c
b++
}
C = a / b
return C
}
# 速度係数 vr を仮定したキャパシタンスの周波数特性
function cap2(vr, M, i, w, c, cm, bl) {
for (i = 0; i < M; i++) {
w = 2 * PI * fm[i]
cm = -xm[i] / w / (rm[i] * rm[i] + xm[i] * xm[i])
bl = w / vr / 2.99792458e8 * l
c = cm / sin(bl) * cos(bl) * bl
X[i] = i * 1e6
Y[i] = c
}
ve = reg0(M)
return ve # 分散
}
# Y[i] (i = 0, 1, .. N-1) の平均と分散
function reg0(N, i, T, Q, S, V) {
T = Q = 0
for (i = 0; i < N; i++) {
T += Y[i]
Q += Y[i] * Y[i]
}
S = Q - T * T / N # 変動
V = S / (N - 1) # 分散
# average = T / N # 平均
return V
}
# 回帰直線を求める
# Y[i] = b0 + b1 * X[i] (i = 1, 2, .. N)
function reg1(N, \
i, x1, y1, x2, y2, tx, ty, qx, qy, sx, sy, cxy, a, b, ve) {
tx = ty = qx = qy = cxy = offset = 0
for (i = 0; i < N; i++) {
tx += X[i]
ty += Y[i]
qx += X[i] * X[i]
qy += Y[i] * Y[i]
cxy += X[i] * Y[i]
}
sx = qx - tx * tx / N
sy = qy - ty * ty / N
cxy -= tx * ty / N
b1 = cxy / sx
b0 = ty / N - a * (tx / N)
ve = (sy - cxy * cxy / sx) / (N - 2)
return ve # 分散
}
# ２変数の重回帰分析
# y = b0 + b1 * X1[] + b2 * X2[]
function reg2(N, T1, T2, Ty, S11, S22, S12, Syy, S1y, S2y) {
T1 = T2 = Ty = S11 = S22 = S12 = Syy = S1y = S2y = 0
for (i = 0; i < N; i++) {
T1 += X1[i]
T2 += X2[i]
Ty += Y[i]
S11 += X1[i] * X1[i]
S22 += X2[i] * X2[i]
S12 += X1[i] * X2[i]
Syy += Y[i] * Y[i]
S1y += X1[i] * Y[i]
S2y += X2[i] * Y[i]
}
Ty /= N; T1 /= N; T2 /= N # 平均
# T1 = T2 = Ty = 0 # y = b1 * X1 + b2 * X2 (b0 = 0) model
S11 -= N * T1 * T1 # 変動
S22 -= N * T2 * T2
S12 -= N * T1 * T2 # 共変動
Syy -= N * Ty * Ty
S1y -= N * T1 * Ty
S2y -= N * T2 * Ty
b1 = (S1y * S22 - S2y * S12) / (S11 * S22 - S12 * S12)
b2 = (S2y * S11 - S1y * S12) / (S11 * S22 - S12 * S12)
b0 = Ty - b1 * T1 - b2 * T2
return
for (i = 0; i < N; i++) {
y = b1 * X1[i] + b2 * X2[i]
y1 = b1 * X1[i]
y0 = y + b0
print Y[i], y, y/y0, y1/y
}
}
# 入力インピーダンス測定
function meas(f, a, b, w) {
w = 2 * PI * f
a = _a * sqrt(f) + _b * f
b = w / 2.99792458e8 * sqrt(2.3) + a
Zoc(_Z0, a, b, l)
}
# キャパシタンス測定
function meas_c(f, a, b, w) {
w = 2 * PI * f
a = _a * sqrt(f) + _b * f
b = w / 2.99792458e8 / _v
Zoc(_Z0, a, b, l)
}
# 次の共振周波数を探す
# 入力インピーダンスの虚部が 0 になる点を探しています
function find_next(f1, f2, val, i, f, w, a, b, y1, y2) {
if (1) {
meas(f1); y1 = y
meas(f2); y2 = y
for (i = 0; i < 50; i++) {
f = (f1 + f2) * 0.5
meas(f)
if (y == 0.0 || (f2 - f1) / f < 1e-6)
break
else if (y1 * y > 0) {
f1 = f; y1 = y
}
else {
f2 = f; y2 = y
}
}
} else {
for (i = 0; i < 50; i++) {
f = (f1 + f2) * 0.5
meas(f)
#print i, x, y, f1, f2
if (y == 0.0 || (f2 - f1) / f < 1e-6)
break
else if (y < 0)
f1 = f
else
f2 = f
}
}
return f
}
# open circuit 1/4 wave length frequency
# 長さ l 速度係数 Vr 近くの終端開放線路の 1/4 波長周波数を見付ける
function find_f0o(l, Vr, i, f, f1, f2, w, a, b)
{
f0 = 2.99792458e8 * Vr / l * 0.25
f1 = f0 * 0.75; f2 = f0 * 1.25
if (1) {
f = find_next(f1, f2)
} else {
for (i = 0; i < 50; i++) {
f = (f1 + f2) * 0.5
meas(f)
if (y == 0.0 || (f2 - f1) / f < 1e-6)
break
else if (y < 0)
f1 = f
else
f2 = f
}
}
return f
}
# short circuit 1/4 wave length frequency
function find_f0s(l, Vr, i, f, f1, f2, w, a, b)
{
f0 = 2.99792458e8 * Vr / l * 0.25
f1 = f0 * 0.75; f2 = f0 * 1.25
for (i = 0; i < 50; i++) {
f = (f1 + f2) * 0.5
w = 2 * PI * f
a = _a * sqrt(f) + _b * f
b = w / (2.99792458e8 * Vr) + a
Zsc(_Z0, a, b, l)
if (y == 0.0 || (f2 - f1) / (f1 + f2) < 1e-12)
break
else if (y > 0)
f1 = f
else
f2 = f
}
return f
}
# Zin = (Zt/Z0 + tanh(γ*l) / (1 + (Zt/Z0)*tanh(γ*l)) * Z0
function Zin(z0, x1, y1, l, x2, y2, x3, y3, x4, y4) {
Tanh(x1 * l, y1 * l)
x3 = x
y3 = y
Mul(x, y, x2 / z0, y2 / z0)
x4 = x + 1
y4 = y
Add(x3, y3, x2 / z0, y2 / z0)
Div(x, y, x4, y4)
x *= z0
y *= z0
}
# Zsc = Z0 * tanh((α + j * β) * l)
function Zsc(z0, x1, y1, l) {
Tanh(x1 * l, y1 * l)
Mul(x, y, z0, 0)
}
# Zoc = Z0 * coth((α + j * β) * l)
function Zoc(z0, x1, y1, l, w) {
Tanh(x1 * l, y1 * l)
Div(z0, 0, x, y)
}
function atanh(x) {
return 0.5 * log((1 + x) / ( 1 - x))
}
# complex functions
function Add(x1, y1, x2, y2) {
x = x1 + x2
y = y1 + y2
}
function Sub(x1, y1, x2, y2) {
x = x1 - x2
y = y1 - y2
}
function Mul(x1, y1, x2, y2) {
x = x1 * x2 - y1 * y2
y = x1 * y2 + y1 * x2
}
function Div(x1, y1, x2, y2, r2) {
r2 = x2 * x2 + y2 * y2
x = (x1 * x2 + y1 * y2) / r2
y = (y1 * x2 - x1 * y2) / r2
}
function Exp(x1, y1, a) {
a = exp(x1)
x = a * cos(y1)
y = a * sin(y1)
}
function Log(x1, y1) {
x = 0.5 * log(x1 * x1 + y1 * y1)
y = atan2(y1, x1)
}
function Sqrt(x1, y1, r, w) {
x = 0.5 * log(x1 * x1 + y1 * y1) * 0.5
y = atan2(y1, x1) * 0.5
a = exp(x)
x = a * cos(y)
y = a * sin(y)
}
function Pow(x1, y1, x2, y2, x3, y3) {
Log(x1, y1)
x3 = x
y3 = y
Mul(x3, y3, x2, y2)
Exp(x, y)
}
function Sin(x1, y1, e, f) {
e = exp(y1)
f = 1 / e
x = 0.5 * cos(x1) * (e - f)
y = 0.5 * sin(x1) * (e + f)
}
function Cos(x1, y1, e, f) {
e = exp(y1)
f = 1 / e
x = 0.5 * cos(x1) * (f + e)
y = 0.5 * sin(x1) * (f - e)
}
function Tan(x1, y1, e, f, d) {
e = exp(2 * y1)
f = 1 / e
d = cos(2 * x1) + 0.5 * (e + f)
x = sin(2 * x1) / d
y = 0.5 * (e - f) / d
}
function Sinh(x1, y1, e, f) {
e = exp(x1)
f = 1 / e
x = 0.5 * (e - f) * cos(y1)
y = 0.5 * (e + f) * sin(y1)
}
function Cosh(x1, y1, e, f) {
e = exp(x1)
f = 1 / e
x = 0.5 * (e + f) * cos(y1)
y = 0.5 * (e - f) * sin(y1)
}
function Tanh(x1, y1, e, f, d) {
e = exp(2 * x1)
f = 1 / e
d = 0.5 * (e + f) + cos(2 * y1)
if (d == 0) {
e = exp(x1)
f = 1 / e
x = 0
y = (e + f) * sin(y1)
}
else {
x = 0.5 * (e - f) / d
y = sin(2 * y1) / d
}
}
function Atanh(x1, y1) {
Div(1 + x1, y1, 1 - x1, -y1)
Log(x, y)
x *= 0.5
y *= 0.5
}
function Abs(x, y) {
return sqrt(x * x + y * y)
}
function Arg(x, y) {
return atan2(y, x)
}
'
高周波インピーダンス測定機は安価な機械ではありませんから、 「そんな高価なものは使えないよ」という状況もたくさんありますが、 二次定数測定には実にたくさんの方法があって、 いろいろな工夫ができます。 私もかつてはディップメータ、Q メータ、 VSWR ブリッジと発信機とか、 発信機とオシロスコープを組み合わせた即席TDRとか、 ケーブル終端に可変抵抗を接続して、 その電圧の出力変化から特性インピーダンスを求めるとか、 いろいろ工夫していました。 こういった工夫と実験もまた、 深い理解を得るためのよい手段になります。 「そうか」「わかった」という理解は感動の源泉です。
平林浩一, 2016-05