PSK（Phase Shift Keying: 位相偏移変調）は、搬送波の位相にディジタル情報を載せる変調方式です。AM変調やFM変調がアナログ変調であったのに対し、PSKはディジタル通信の中核をなす変調方式であり、Wi-Fi、LTE、衛星通信、深宇宙通信など、現代の通信システムで広く使用されています。
PSKの大きな利点は、振幅を一定に保てるため非線形増幅器でも歪みなく使用できる点です。本記事では、BPSK、QPSK、8PSKの理論とビット誤り率（BER）を数学的に導出し、Pythonでシミュレーションを行います。
本記事の内容
- PSKの基本原理（位相にビット情報を載せる仕組み）
- BPSK: 信号点配置とBER導出
- QPSK: 直交成分分解とグレイ符号化
- 8PSK: 帯域効率と誤り率のトレードオフ
- 差動符号化（DPSK）
- Pythonでの星座図・BER理論曲線・雑音シミュレーション
前提知識
この記事を読む前に、以下の記事を読んでおくと理解が深まります。
PSKの基本原理
ディジタル変調の考え方
ディジタル変調では、ディジタルデータ（0と1のビット列）を搬送波のパラメータの変化として表現します。PSKでは位相を変化させます。
$M$ 値PSK（M-PSK）では、$\log_2 M$ ビットのデータを1つのシンボルに割り当て、$M$ 個の位相状態のいずれかに対応させます。
$$ s_i(t) = A_c \cos\!\left(2\pi f_c t + \phi_i\right), \quad i = 0, 1, \ldots, M-1 $$
ここで $\phi_i$ は $i$ 番目のシンボルに割り当てられた位相です。等間隔に位相を配置する場合、
$$ \phi_i = \frac{2\pi i}{M}, \quad i = 0, 1, \ldots, M-1 $$
IQ平面表現
PSK信号は、同相成分（In-phase: I）と直交成分（Quadrature: Q）に分解して表現するのが便利です。
$$ s_i(t) = A_c \cos\phi_i \cos(2\pi f_c t) – A_c \sin\phi_i \sin(2\pi f_c t) $$
$$ = I_i \cos(2\pi f_c t) – Q_i \sin(2\pi f_c t) $$
ここで
$$ I_i = A_c \cos\phi_i, \quad Q_i = A_c \sin\phi_i $$
各シンボルは、IQ平面（複素平面）上の点 $I_i + jQ_i = A_c e^{j\phi_i}$ として表現されます。この図を星座図（constellation diagram）と呼びます。
M-PSKの信号点は、IQ平面上の半径 $A_c$ の円周上に等間隔に配置されます。
BPSK（Binary PSK）
信号点配置
BPSK（$M = 2$）では、1シンボルに1ビットを割り当てます。位相は2つです。
$$ \phi_0 = 0 \quad (\text{ビット 0}), \quad \phi_1 = \pi \quad (\text{ビット 1}) $$
信号は次のようになります。
$$ s_0(t) = A_c \cos(2\pi f_c t) \quad (\text{ビット 0}) $$
$$ s_1(t) = A_c \cos(2\pi f_c t + \pi) = -A_c \cos(2\pi f_c t) \quad (\text{ビット 1}) $$
IQ平面では、2つの信号点は実軸上の $(A_c, 0)$ と $(-A_c, 0)$ に位置します。2点間の距離は $2A_c$ です。
BPSK の BER 導出
シンボルエネルギーとビットエネルギー
1シンボルの持続時間を $T_b$（ビット周期 = シンボル周期）とすると、1ビットあたりのエネルギーは
$$ E_b = \frac{A_c^2}{2} T_b $$
ここで $A_c^2/2$ は $s(t) = \pm A_c \cos(2\pi f_c t)$ の平均電力です。
判定規則
受信信号を $r(t) = s_i(t) + n(t)$ とします。$n(t)$ は平均0、電力スペクトル密度 $N_0/2$ の加法性白色ガウス雑音（AWGN）です。
整合フィルタ（matched filter）の出力は、
$$ z = \int_0^{T_b} r(t) \cdot \sqrt{\frac{2}{T_b}} \cos(2\pi f_c t) \, dt $$
基底関数を $\phi_1(t) = \sqrt{2/T_b} \cos(2\pi f_c t)$ とすると、
ビット0が送信された場合: $s_0(t) = \sqrt{E_b} \cdot \phi_1(t)$ より $z_s = \sqrt{E_b}$
ビット1が送信された場合: $s_1(t) = -\sqrt{E_b} \cdot \phi_1(t)$ より $z_s = -\sqrt{E_b}$
雑音成分 $z_n$ は平均0、分散 $N_0/2$ のガウス分布に従います。
$$ z = z_s + z_n, \quad z_n \sim \mathcal{N}(0, N_0/2) $$
判定は $z \gtrless 0$ で行います。
BER の計算
ビット0が送信されたとき（$z_s = \sqrt{E_b}$）に誤る確率は、
$$ P_e = P(z < 0 \mid s_0) = P\!\left(z_n < -\sqrt{E_b}\right) $$
$z_n \sim \mathcal{N}(0, N_0/2)$ を標準正規分布に変換します。$u = z_n / \sqrt{N_0/2}$ とおくと $u \sim \mathcal{N}(0, 1)$ であり、
$$ P_e = P\!\left(u < \frac{-\sqrt{E_b}}{\sqrt{N_0/2}}\right) = P\!\left(u < -\sqrt{\frac{2E_b}{N_0}}\right) $$
Q関数 $Q(x) = P(u > x)$ を用いて、
$$ P_e = Q\!\left(\sqrt{\frac{2E_b}{N_0}}\right) $$
対称性より、ビット1が送信された場合も同じ誤り確率となります。したがって、
$$ \boxed{ P_b^{\text{BPSK}} = Q\!\left(\sqrt{\frac{2E_b}{N_0}}\right) } $$
ここで Q関数は次のように定義されます。
$$ Q(x) = \frac{1}{\sqrt{2\pi}} \int_x^{\infty} e^{-t^2/2} \, dt = \frac{1}{2} \text{erfc}\!\left(\frac{x}{\sqrt{2}}\right) $$
QPSK（Quadrature PSK）
信号点配置
QPSK（$M = 4$）では、2ビットを1シンボルに割り当てます。4つの位相は
$$ \phi_i = \frac{\pi}{4} + \frac{\pi i}{2}, \quad i = 0, 1, 2, 3 $$
IQ平面上の信号点は
$$ s_i = A_c e^{j\phi_i} $$
具体的に：
| シンボル $i$ | ビット（グレイ） | 位相 $\phi_i$ | $(I_i, Q_i)$ |
|---|---|---|---|
| 0 | 00 | $\pi/4$ | $(A_c/\sqrt{2}, A_c/\sqrt{2})$ |
| 1 | 01 | $3\pi/4$ | $(-A_c/\sqrt{2}, A_c/\sqrt{2})$ |
| 2 | 11 | $5\pi/4$ | $(-A_c/\sqrt{2}, -A_c/\sqrt{2})$ |
| 3 | 10 | $7\pi/4$ | $(A_c/\sqrt{2}, -A_c/\sqrt{2})$ |
直交成分分解
QPSK信号は、I成分とQ成分の2つの独立なBPSK信号の和として解釈できます。
$$ s(t) = I_i \cos(2\pi f_c t) – Q_i \sin(2\pi f_c t) $$
$I_i = \pm A_c / \sqrt{2}$, $Q_i = \pm A_c / \sqrt{2}$ ですから、I軸とQ軸でそれぞれ独立にBPSK判定を行えます。
グレイ符号化
隣接するシンボル間で1ビットしか異ならないようにビットを割り当てる方式をグレイ符号化（Gray coding）と呼びます。
雑音によってシンボル判定を誤る場合、最も確率が高いのは隣接シンボルへの誤りです。グレイ符号化により、隣接シンボルへの誤りは1ビットの誤りに対応するため、ビット誤り率が最小化されます。
QPSK の BER 導出
QPSKをI成分とQ成分の2つの独立なBPSKとみなします。
各軸上の信号点距離は $2A_c/\sqrt{2} = \sqrt{2}A_c$ です。各軸のビットエネルギーを求めます。
1シンボルのエネルギーは $E_s = A_c^2 T_s / 2$（$T_s$ はシンボル周期）です。QPSKでは1シンボルに2ビット含まれるので、
$$ E_b = \frac{E_s}{2} = \frac{A_c^2 T_s}{4} $$
各軸（I, Q）の信号は $\pm A_c/\sqrt{2}$ のBPSKです。各軸の平均エネルギーは
$$ E_{\text{axis}} = \frac{(A_c/\sqrt{2})^2}{2} T_s = \frac{A_c^2 T_s}{4} = E_b $$
したがって、各軸のBER は
$$ P_{b,\text{axis}} = Q\!\left(\sqrt{\frac{2E_{\text{axis}}}{N_0}}\right) = Q\!\left(\sqrt{\frac{2E_b}{N_0}}\right) $$
グレイ符号化の下では、全体のBERは各軸のBERの平均です。I軸とQ軸は独立なので、
$$ \boxed{ P_b^{\text{QPSK}} = Q\!\left(\sqrt{\frac{2E_b}{N_0}}\right) } $$
驚くべきことに、QPSKのBERはBPSKと全く同じです。QPSKは帯域幅あたりの伝送効率（帯域効率）がBPSKの2倍でありながら、ビット誤り率は同等です。これがQPSKが広く使われる理由です。
帯域効率の比較
ビットレートを $R_b$、帯域幅を $B$ とすると、
| 方式 | シンボルレート $R_s$ | 帯域幅 $B$ | 帯域効率 $R_b/B$ |
|---|---|---|---|
| BPSK | $R_b$ | $2R_b$ | 0.5 bit/s/Hz |
| QPSK | $R_b/2$ | $R_b$ | 1 bit/s/Hz |
QPSKはシンボルレートがBPSKの半分なので占有帯域幅も半分になり、帯域効率が2倍に改善されます。
8PSK
信号点配置
8PSK（$M = 8$）では、3ビットを1シンボルに割り当てます。8つの位相は
$$ \phi_i = \frac{2\pi i}{8} = \frac{\pi i}{4}, \quad i = 0, 1, \ldots, 7 $$
帯域効率は $\log_2 8 = 3$ bit/s/Hz です。
8PSK の BER 近似
M-PSKの一般的なBER近似式（グレイ符号化、$M \geq 4$）は、
$$ P_b \approx \frac{2}{\log_2 M} Q\!\left(\sqrt{2 \log_2 M \cdot \frac{E_b}{N_0}} \sin\frac{\pi}{M}\right) $$
この式は、隣接シンボルへの誤り確率がドミナントであるという近似に基づきます。
導出の概要
隣接する2シンボル間の距離は
$$ d_{\min} = 2A_c \sin\!\left(\frac{\pi}{M}\right) $$
M-PSKにおけるシンボル誤り率（SER）は、隣接シンボルへの誤り確率の2倍で近似できます（各シンボルは2つの隣接シンボルを持つため）。
$$ P_s \approx 2Q\!\left(\frac{d_{\min}/2}{\sqrt{N_0/2}}\right) = 2Q\!\left(\frac{A_c \sin(\pi/M)}{\sqrt{N_0/2}}\right) $$
シンボルエネルギー $E_s = A_c^2 T_s / 2$ とビットエネルギー $E_b = E_s / \log_2 M$ の関係から $A_c^2 = 2E_s / T_s = 2\log_2 M \cdot E_b / T_s$ を代入し整理すると、
$$ P_s \approx 2Q\!\left(\sqrt{2\log_2 M \cdot \frac{E_b}{N_0}} \sin\frac{\pi}{M}\right) $$
グレイ符号化の下では、1シンボル誤りが1ビット誤りに対応するので、
$$ P_b \approx \frac{P_s}{\log_2 M} = \frac{2}{\log_2 M} Q\!\left(\sqrt{2\log_2 M \cdot \frac{E_b}{N_0}} \sin\frac{\pi}{M}\right) $$
8PSKの場合（$M = 8$）、
$$ P_b^{\text{8PSK}} \approx \frac{2}{3} Q\!\left(\sqrt{6 \cdot \frac{E_b}{N_0}} \sin\frac{\pi}{8}\right) $$
$\sin(\pi/8) \approx 0.3827$ です。
帯域効率と誤り率のトレードオフ
$M$ を増やすと帯域効率は $\log_2 M$ に比例して向上しますが、信号点間距離が縮小するため、同じ $E_b/N_0$ に対するBERが悪化します。
| 方式 | 帯域効率 [bit/s/Hz] | 同BERに必要な$E_b/N_0$の増加量 |
|---|---|---|
| BPSK | 0.5 | 基準 |
| QPSK | 1.0 | 0 dB |
| 8PSK | 1.5 | 約 3.5 dB |
| 16PSK | 2.0 | 約 8.2 dB |
帯域効率と雑音耐性は本質的にトレードオフの関係にあります。
差動符号化（DPSK）
位相不確定性の問題
コヒーレント検波では、受信側で搬送波の絶対位相を知る必要があります。しかし、キャリア再生回路は $\pi$ の不確定性（BPSKの場合）や $\pi/2$ の不確定性（QPSKの場合）を持つ場合があります。
DPSK（Differential PSK）
DPSKでは、情報をシンボルの絶対位相ではなく、連続するシンボル間の位相差に載せます。
DBPSKの場合、送信する位相差は
$$ \Delta\phi_k = \phi_k – \phi_{k-1} = \begin{cases} 0 & (\text{ビット 0}) \\ \pi & (\text{ビット 1}) \end{cases} $$
復調は、現在のシンボルと前のシンボルの位相差を検出するだけなので、絶対位相の推定が不要です。
DBPSKのBER
$$ P_b^{\text{DBPSK}} = \frac{1}{2} e^{-E_b/N_0} $$
BPSKに比べて若干BERが悪化しますが（約1 dBの劣化）、位相不確定性の問題を回避できます。
DQPSKのBER近似
$$ P_b^{\text{DQPSK}} \approx Q\!\left(\sqrt{\frac{2E_b}{N_0}} \cdot 2\sin\frac{\pi}{4\sqrt{2}}\right) $$
これも通常のQPSKに比べて若干の劣化がありますが、実装が簡単という利点があります。
Pythonでの実装
各種PSKの星座図
import numpy as np
import matplotlib.pyplot as plt
def generate_constellation(M):
"""M-PSKの星座図を生成する"""
phases = np.array([2 * np.pi * i / M for i in range(M)])
symbols = np.exp(1j * phases)
return symbols, phases
fig, axes = plt.subplots(2, 2, figsize=(12, 12))
axes = axes.flatten()
M_list = [2, 4, 8, 16]
titles = ['BPSK (M=2)', 'QPSK (M=4)', '8PSK (M=8)', '16PSK (M=16)']
for idx, (M, title) in enumerate(zip(M_list, titles)):
symbols, phases = generate_constellation(M)
ax = axes[idx]
# 単位円
theta_circle = np.linspace(0, 2*np.pi, 200)
ax.plot(np.cos(theta_circle), np.sin(theta_circle), 'k--',
linewidth=0.5, alpha=0.3)
# 信号点
ax.scatter(symbols.real, symbols.imag, s=100, c='red',
zorder=5, edgecolors='darkred', linewidths=1.5)
# ビットラベル（グレイ符号）
bits_per_symbol = int(np.log2(M))
for i, sym in enumerate(symbols):
# グレイ符号化
gray = i ^ (i >> 1)
label = format(gray, f'0{bits_per_symbol}b')
ax.annotate(label, (sym.real, sym.imag),
textcoords="offset points",
xytext=(10, 10), fontsize=9,
color='blue', fontweight='bold')
# 判定領域の境界線
for i in range(M):
boundary_angle = phases[i] + np.pi / M
ax.plot([0, 1.3*np.cos(boundary_angle)],
[0, 1.3*np.sin(boundary_angle)],
'g--', linewidth=0.8, alpha=0.5)
ax.set_xlim(-1.5, 1.5)
ax.set_ylim(-1.5, 1.5)
ax.set_aspect('equal')
ax.grid(True, alpha=0.3)
ax.axhline(y=0, color='k', linewidth=0.5)
ax.axvline(x=0, color='k', linewidth=0.5)
ax.set_xlabel('In-Phase (I)')
ax.set_ylabel('Quadrature (Q)')
ax.set_title(title)
plt.suptitle('PSK Constellation Diagrams with Gray Coding', fontsize=14)
plt.tight_layout()
plt.show()
星座図上の信号点が単位円上に等間隔に配置され、グレイ符号化により隣接シンボルが1ビットだけ異なることが確認できます。
BER理論曲線
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import erfc
def Q_func(x):
"""Q関数の計算"""
return 0.5 * erfc(x / np.sqrt(2))
# Eb/N0の範囲
EbN0_dB = np.linspace(0, 16, 200)
EbN0 = 10**(EbN0_dB / 10)
# BPSK (= QPSK) の理論BER
ber_bpsk = Q_func(np.sqrt(2 * EbN0))
# 8PSK の理論BER（近似式）
ber_8psk = (2/3) * Q_func(np.sqrt(6 * EbN0) * np.sin(np.pi/8))
# 16PSK の理論BER（近似式）
ber_16psk = (2/4) * Q_func(np.sqrt(8 * EbN0) * np.sin(np.pi/16))
# DBPSK の理論BER
ber_dbpsk = 0.5 * np.exp(-EbN0)
fig, ax = plt.subplots(figsize=(10, 8))
ax.semilogy(EbN0_dB, ber_bpsk, 'b-', linewidth=2, label='BPSK / QPSK')
ax.semilogy(EbN0_dB, ber_8psk, 'r-', linewidth=2, label='8PSK')
ax.semilogy(EbN0_dB, ber_16psk, 'g-', linewidth=2, label='16PSK')
ax.semilogy(EbN0_dB, ber_dbpsk, 'b--', linewidth=2, label='DBPSK')
ax.set_xlabel('$E_b/N_0$ [dB]', fontsize=12)
ax.set_ylabel('Bit Error Rate (BER)', fontsize=12)
ax.set_title('BER vs $E_b/N_0$ for M-PSK', fontsize=14)
ax.legend(fontsize=11)
ax.grid(True, which='both', alpha=0.3)
ax.set_xlim(0, 16)
ax.set_ylim(1e-6, 1)
# BER = 10^-3 のラインを引く
ax.axhline(y=1e-3, color='gray', linestyle=':', alpha=0.5)
ax.annotate('BER = $10^{-3}$', xy=(1, 1e-3), fontsize=10, color='gray')
plt.tight_layout()
plt.show()
# BER = 10^-3 に必要な Eb/N0 を表示
print("=== BER = 10^-3 に必要な Eb/N0 ===")
target_ber = 1e-3
for name, ber in [('BPSK/QPSK', ber_bpsk), ('8PSK', ber_8psk),
('16PSK', ber_16psk)]:
idx = np.argmin(np.abs(ber - target_ber))
print(f" {name:12s}: Eb/N0 = {EbN0_dB[idx]:.1f} dB")
QPSKとBPSKのBER曲線が完全に一致すること、そして $M$ が増えるにつれて同じBERを達成するために必要な $E_b/N_0$ が増加することが確認できます。
雑音付加シミュレーション（モンテカルロ法）
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import erfc
def Q_func(x):
"""Q関数"""
return 0.5 * erfc(x / np.sqrt(2))
def psk_simulation(M, EbN0_dB_range, num_symbols=100000):
"""
M-PSKのモンテカルロBERシミュレーション
"""
bits_per_symbol = int(np.log2(M))
ber_sim = []
# 星座図の生成
constellation = np.exp(1j * 2 * np.pi * np.arange(M) / M)
# グレイ符号化のマッピング
gray_map = np.zeros(M, dtype=int)
for i in range(M):
gray_map[i] = i ^ (i >> 1)
# 逆マッピング（グレイ符号 -> 自然数）
inv_gray = np.zeros(M, dtype=int)
for i in range(M):
inv_gray[gray_map[i]] = i
for EbN0_dB in EbN0_dB_range:
EbN0 = 10**(EbN0_dB / 10)
Es = EbN0 * bits_per_symbol # シンボルエネルギー/N0
# ランダムなシンボルインデックスを生成
sym_idx = np.random.randint(0, M, num_symbols)
# 送信信号（グレイ符号化された星座点）
tx_symbols = constellation[gray_map[sym_idx]]
# AWGN雑音の付加
# シンボルエネルギーをEsに正規化（|constellation|=1なのでそのまま）
noise_std = np.sqrt(1 / (2 * Es))
noise = noise_std * (np.random.randn(num_symbols)
+ 1j * np.random.randn(num_symbols))
rx_symbols = tx_symbols + noise
# 最尤判定（最近傍シンボルを選択）
distances = np.abs(rx_symbols[:, None] - constellation[gray_map][None, :])
detected_idx = np.argmin(distances, axis=1)
# ビット誤り数のカウント
tx_bits = np.array([list(format(gray_map[s], f'0{bits_per_symbol}b'))
for s in sym_idx], dtype=int)
rx_bits = np.array([list(format(gray_map[d], f'0{bits_per_symbol}b'))
for d in detected_idx], dtype=int)
bit_errors = np.sum(tx_bits != rx_bits)
total_bits = num_symbols * bits_per_symbol
ber_sim.append(bit_errors / total_bits)
return np.array(ber_sim)
# シミュレーション実行
EbN0_dB_range = np.arange(0, 17, 1)
EbN0_theory = np.linspace(0, 16, 200)
EbN0_lin = 10**(EbN0_theory / 10)
print("BPSKシミュレーション実行中...")
ber_bpsk_sim = psk_simulation(2, EbN0_dB_range)
print("QPSKシミュレーション実行中...")
ber_qpsk_sim = psk_simulation(4, EbN0_dB_range)
print("8PSKシミュレーション実行中...")
ber_8psk_sim = psk_simulation(8, EbN0_dB_range)
# 理論値
ber_bpsk_theory = Q_func(np.sqrt(2 * EbN0_lin))
ber_8psk_theory = (2/3) * Q_func(np.sqrt(6 * EbN0_lin) * np.sin(np.pi/8))
# プロット
fig, ax = plt.subplots(figsize=(10, 8))
# 理論曲線
ax.semilogy(EbN0_theory, ber_bpsk_theory, 'b-', linewidth=2,
label='BPSK/QPSK (Theory)')
ax.semilogy(EbN0_theory, ber_8psk_theory, 'r-', linewidth=2,
label='8PSK (Theory)')
# シミュレーション結果
mask_bpsk = ber_bpsk_sim > 0
mask_qpsk = ber_qpsk_sim > 0
mask_8psk = ber_8psk_sim > 0
ax.semilogy(EbN0_dB_range[mask_bpsk], ber_bpsk_sim[mask_bpsk], 'bs',
markersize=8, label='BPSK (Sim)')
ax.semilogy(EbN0_dB_range[mask_qpsk], ber_qpsk_sim[mask_qpsk], 'g^',
markersize=8, label='QPSK (Sim)')
ax.semilogy(EbN0_dB_range[mask_8psk], ber_8psk_sim[mask_8psk], 'ro',
markersize=8, label='8PSK (Sim)')
ax.set_xlabel('$E_b/N_0$ [dB]', fontsize=12)
ax.set_ylabel('Bit Error Rate (BER)', fontsize=12)
ax.set_title('PSK BER: Theory vs Monte Carlo Simulation', fontsize=14)
ax.legend(fontsize=10)
ax.grid(True, which='both', alpha=0.3)
ax.set_xlim(0, 16)
ax.set_ylim(1e-5, 1)
plt.tight_layout()
plt.show()
シミュレーション結果が理論曲線とよく一致していることが確認できます。BPSKとQPSKのBERが等しいこと、8PSKが同じ $E_b/N_0$ に対してBERが悪化することも実証されています。
雑音環境下の星座図
import numpy as np
import matplotlib.pyplot as plt
def plot_noisy_constellation(M, EbN0_dB, num_symbols=3000, ax=None):
"""雑音環境下の星座図をプロットする"""
bits_per_symbol = int(np.log2(M))
EbN0 = 10**(EbN0_dB / 10)
Es = EbN0 * bits_per_symbol
# 星座図
constellation = np.exp(1j * 2 * np.pi * np.arange(M) / M)
gray_map = np.array([i ^ (i >> 1) for i in range(M)])
# ランダムシンボル生成
sym_idx = np.random.randint(0, M, num_symbols)
tx = constellation[gray_map[sym_idx]]
# AWGN付加
noise_std = np.sqrt(1 / (2 * Es))
noise = noise_std * (np.random.randn(num_symbols)
+ 1j * np.random.randn(num_symbols))
rx = tx + noise
if ax is None:
fig_local, ax = plt.subplots(figsize=(6, 6))
ax.scatter(rx.real, rx.imag, s=1, alpha=0.3, c='blue')
ax.scatter(constellation.real, constellation.imag,
s=100, c='red', zorder=5, edgecolors='darkred')
# 単位円
theta_c = np.linspace(0, 2*np.pi, 200)
ax.plot(np.cos(theta_c), np.sin(theta_c), 'k--', linewidth=0.5, alpha=0.3)
ax.set_xlim(-2, 2)
ax.set_ylim(-2, 2)
ax.set_aspect('equal')
ax.grid(True, alpha=0.3)
ax.set_xlabel('In-Phase (I)')
ax.set_ylabel('Quadrature (Q)')
fig, axes = plt.subplots(3, 3, figsize=(15, 15))
configs = [
(4, 3), (4, 7), (4, 12),
(8, 3), (8, 7), (8, 12),
(16, 3), (16, 7), (16, 12),
]
for idx, (M, snr) in enumerate(configs):
row, col = divmod(idx, 3)
plot_noisy_constellation(M, snr, ax=axes[row, col])
name = {2:'BPSK', 4:'QPSK', 8:'8PSK', 16:'16PSK'}[M]
axes[row, col].set_title(f'{name}, $E_b/N_0$ = {snr} dB')
plt.suptitle('Noisy Constellation Diagrams', fontsize=14, y=1.01)
plt.tight_layout()
plt.show()
$E_b/N_0$ が低い（雑音が大きい）環境では、信号点が大きく拡散し、隣接シンボルとの判定境界を越える確率が増加することが視覚的に理解できます。$M$ が大きいほどシンボル間距離が近く、雑音の影響を受けやすいことも確認できます。
まとめ
本記事では、PSK（位相偏移変調）の理論と星座図について解説しました。
- M-PSKは搬送波の位相に $\log_2 M$ ビットの情報を載せるディジタル変調方式である
- BPSKのBERは $P_b = Q(\sqrt{2E_b/N_0})$ であり、AWGN環境での最適性を持つ
- QPSKはBPSKと同じBERでありながら、帯域効率が2倍（1 bit/s/Hz）に改善される
- 8PSK以上では帯域効率は向上するが、同じBERに必要な $E_b/N_0$ が増加するトレードオフがある
- グレイ符号化により、隣接シンボルへの誤りが1ビット誤りに抑えられ、BERが最小化される
- DPSK（差動符号化）は位相不確定性を回避できるが、約1 dBのBER劣化がある
次のステップとして、以下の記事も参考にしてください。