スペクトル拡散通信 — DS-SSとFH-SSの統一的理解

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)

# パラメータ
N_bits = 8             # 情報ビット数
N_chips_per_bit = 31   # 1ビットあたりのチップ数（処理利得）
samples_per_chip = 4   # 1チップあたりのサンプル数

# 情報ビット生成（BPSK: +1/-1）
data_bits = 2 * np.random.randint(0, 2, N_bits) - 1

# PN符号生成（M系列の代わりにランダム）
pn_code = 2 * np.random.randint(0, 2, N_chips_per_bit) - 1

# 拡散: 各ビットをPN符号で拡散
spread_signal = np.array([])
for bit in data_bits:
    spread_signal = np.append(spread_signal, bit * pn_code)

# サンプリング（波形生成用）
def upsample(signal, factor):
    return np.repeat(signal, factor)

data_waveform = upsample(np.repeat(data_bits, N_chips_per_bit), samples_per_chip)
spread_waveform = upsample(spread_signal, samples_per_chip)
pn_waveform = upsample(np.tile(pn_code, N_bits), samples_per_chip)

# 雑音付加（低SNR環境）
snr_db = -15  # 拡散帯域でのSNR
noise_power = np.mean(spread_waveform**2) * 10**(-snr_db/10)
noise = np.sqrt(noise_power) * np.random.randn(len(spread_waveform))
received = spread_waveform + noise

# 逆拡散
despread = received * pn_waveform

# ビット判定（チップ期間で積分）
chip_samples = samples_per_chip
bit_samples = N_chips_per_bit * chip_samples
detected_bits = []
for i in range(N_bits):
    start = i * bit_samples
    end = (i + 1) * bit_samples
    integral = np.sum(despread[start:end])
    detected_bits.append(1 if integral > 0 else -1)
detected_bits = np.array(detected_bits)

# 時間軸
t = np.arange(len(data_waveform))

fig, axes = plt.subplots(5, 1, figsize=(14, 12), sharex=True)

axes[0].plot(t, data_waveform, 'b-', linewidth=1.5)
axes[0].set_ylabel('Amplitude')
axes[0].set_title(f'Data bits: {data_bits}')
axes[0].set_ylim([-1.5, 1.5])
axes[0].grid(True, alpha=0.3)

axes[1].plot(t, pn_waveform, 'g-', linewidth=0.5)
axes[1].set_ylabel('Amplitude')
axes[1].set_title(f'PN code ({N_chips_per_bit} chips/bit)')
axes[1].set_ylim([-1.5, 1.5])
axes[1].grid(True, alpha=0.3)

axes[2].plot(t, spread_waveform, 'm-', linewidth=0.5)
axes[2].set_ylabel('Amplitude')
axes[2].set_title('Spread signal (data x PN)')
axes[2].set_ylim([-1.5, 1.5])
axes[2].grid(True, alpha=0.3)

axes[3].plot(t, received, 'r-', linewidth=0.3, alpha=0.7)
axes[3].set_ylabel('Amplitude')
axes[3].set_title(f'Received signal (SNR_spread = {snr_db} dB)')
axes[3].grid(True, alpha=0.3)

axes[4].plot(t, despread, 'c-', linewidth=0.3, alpha=0.7)
# ビット判定結果を重ねて表示
for i in range(N_bits):
    start = i * bit_samples
    end = (i + 1) * bit_samples
    color = 'blue' if detected_bits[i] == data_bits[i] else 'red'
    axes[4].axhspan(-0.3, 0.3, xmin=start/len(t), xmax=end/len(t), alpha=0.2, color=color)
axes[4].set_ylabel('Amplitude')
axes[4].set_title(f'Despread signal → Detected: {detected_bits} (BER = {np.mean(detected_bits != data_bits):.2f})')
axes[4].set_xlabel('Sample')
axes[4].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('dsss_demo.png', dpi=150, bbox_inches='tight')
plt.show()

print(f"送信ビット: {data_bits}")
print(f"検出ビット: {detected_bits}")
print(f"処理利得: {10*np.log10(N_chips_per_bit):.1f} dB")
print(f"拡散帯域SNR: {snr_db} dB")
print(f"逆拡散後の推定SNR: {snr_db + 10*np.log10(N_chips_per_bit):.1f} dB")

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)

# パラメータ
N_bits = 64
N_chips = 63  # 処理利得
samples_per_chip = 8
fs = 1.0  # 正規化サンプリング周波数

# データ生成と拡散
data = 2 * np.random.randint(0, 2, N_bits) - 1
pn = 2 * np.random.randint(0, 2, N_chips) - 1
spread = np.array([])
for bit in data:
    spread = np.append(spread, bit * pn)

# 波形のアップサンプリング
data_wave = np.repeat(np.repeat(data, N_chips), samples_per_chip)
spread_wave = np.repeat(spread, samples_per_chip)
pn_wave = np.repeat(np.tile(pn, N_bits), samples_per_chip)

# 雑音
noise = 3.0 * np.random.randn(len(spread_wave))
received = spread_wave + noise

# 逆拡散
despread = received * pn_wave

# パワースペクトル密度（PSD）
def compute_psd(signal, nfft=4096):
    f = np.fft.rfftfreq(nfft)
    psd = np.abs(np.fft.rfft(signal, n=nfft))**2 / nfft
    # 平滑化
    window = 5
    psd_smooth = np.convolve(psd, np.ones(window)/window, mode='same')
    return f, 10 * np.log10(psd_smooth + 1e-15)

fig, axes = plt.subplots(2, 2, figsize=(14, 8))

f, psd = compute_psd(data_wave)
axes[0, 0].plot(f, psd, 'b-', linewidth=1)
axes[0, 0].set_title('PSD: Original data')
axes[0, 0].set_ylabel('PSD [dB]')
axes[0, 0].set_xlim([0, 0.5])
axes[0, 0].grid(True, alpha=0.3)

f, psd = compute_psd(spread_wave)
axes[0, 1].plot(f, psd, 'g-', linewidth=1)
axes[0, 1].set_title('PSD: Spread signal')
axes[0, 1].set_ylabel('PSD [dB]')
axes[0, 1].set_xlim([0, 0.5])
axes[0, 1].grid(True, alpha=0.3)

f, psd = compute_psd(received)
axes[1, 0].plot(f, psd, 'r-', linewidth=1)
axes[1, 0].set_title('PSD: Received (spread + noise)')
axes[1, 0].set_xlabel('Normalized frequency')
axes[1, 0].set_ylabel('PSD [dB]')
axes[1, 0].set_xlim([0, 0.5])
axes[1, 0].grid(True, alpha=0.3)

f, psd = compute_psd(despread)
axes[1, 1].plot(f, psd, 'c-', linewidth=1)
axes[1, 1].set_title('PSD: After despreading')
axes[1, 1].set_xlabel('Normalized frequency')
axes[1, 1].set_ylabel('PSD [dB]')
axes[1, 1].set_xlim([0, 0.5])
axes[1, 1].grid(True, alpha=0.3)

plt.tight_layout()
plt.savefig('dsss_spectrum.png', dpi=150, bbox_inches='tight')
plt.show()

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(123)

# パラメータ
fs = 50000           # サンプリング周波数 [Hz]
N_symbols = 8        # シンボル数
T_hop = 0.005        # ホッピング周期 [s]（= シンボル周期: slow FH）
M = 8                # 周波数スロット数
f_base = 2000        # ベース周波数 [Hz]
f_spacing = 500      # スロット間隔 [Hz]
hop_frequencies = f_base + np.arange(M) * f_spacing  # 利用可能周波数

samples_per_hop = int(T_hop * fs)

# データシンボル（BFSK簡易版: +1/-1）
data = 2 * np.random.randint(0, 2, N_symbols) - 1

# ホッピングパターン（PN符号で決定）
hop_pattern = np.random.permutation(M)[:N_symbols] % M

# FH-SS信号生成
fh_signal = np.array([])
time_axis = np.array([])
hop_freq_over_time = np.array([])

for i in range(N_symbols):
    t_hop = np.arange(samples_per_hop) / fs
    f_hop = hop_frequencies[hop_pattern[i]]
    # データで周波数を微小シフト（簡易BFSK）
    f_shift = data[i] * 100  # +/-100 Hz
    segment = np.cos(2 * np.pi * (f_hop + f_shift) * t_hop)
    fh_signal = np.append(fh_signal, segment)
    time_axis = np.append(time_axis, t_hop + i * T_hop)
    hop_freq_over_time = np.append(hop_freq_over_time, np.full(samples_per_hop, f_hop))

# スペクトログラム計算
fig, axes = plt.subplots(3, 1, figsize=(14, 10))

# ホッピングパターン
for i in range(N_symbols):
    f_hop = hop_frequencies[hop_pattern[i]]
    t_start = i * T_hop
    t_end = (i + 1) * T_hop
    axes[0].barh(f_hop/1000, T_hop*1000, left=t_start*1000, height=f_spacing*0.8/1000,
                color=plt.cm.Set3(hop_pattern[i]/M), edgecolor='black', linewidth=0.5)
    axes[0].text((t_start + T_hop/2)*1000, f_hop/1000,
                f'd={data[i]:+d}', ha='center', va='center', fontsize=8)

axes[0].set_ylabel('Frequency [kHz]')
axes[0].set_title('FH-SS Hopping Pattern')
axes[0].set_ylim([f_base/1000 - 0.3, (f_base + M*f_spacing)/1000 + 0.3])
axes[0].grid(True, alpha=0.3)

# 時間波形
axes[1].plot(time_axis * 1000, fh_signal, 'b-', linewidth=0.3)
axes[1].set_ylabel('Amplitude')
axes[1].set_title('FH-SS Time Domain Signal')
axes[1].grid(True, alpha=0.3)

# スペクトログラム
nfft = 256
axes[2].specgram(fh_signal, NFFT=nfft, Fs=fs, noverlap=nfft//2, cmap='inferno')
axes[2].set_ylabel('Frequency [Hz]')
axes[2].set_xlabel('Time [s]')
axes[2].set_title('Spectrogram of FH-SS Signal')
axes[2].set_ylim([f_base - 500, f_base + M * f_spacing + 500])

plt.tight_layout()
plt.savefig('fhss_demo.png', dpi=150, bbox_inches='tight')
plt.show()

print(f"データシンボル: {data}")
print(f"ホッピングパターン: {hop_pattern}")
print(f"使用周波数 [Hz]: {[hop_frequencies[h] for h in hop_pattern]}")

import numpy as np
import matplotlib.pyplot as plt

np.random.seed(42)

# 共通パラメータ
N_bits = 1000
N_chips = 31  # DS-SSの処理利得

# === DS-SS ===
data_dsss = 2 * np.random.randint(0, 2, N_bits) - 1
pn = 2 * np.random.randint(0, 2, N_chips) - 1

# 拡散
spread_chips = np.array([])
for bit in data_dsss:
    spread_chips = np.append(spread_chips, bit * pn)

# 狭帯域ジャミング（単一トーン）
jammer_freq = 0.2  # 正規化周波数
t_chips = np.arange(len(spread_chips))
jammer_power_range = np.arange(0, 30, 2)  # J/S [dB]

ber_dsss = []
for js_db in jammer_power_range:
    js_lin = 10**(js_db/10)
    jammer = np.sqrt(2 * js_lin) * np.cos(2 * np.pi * jammer_freq * t_chips)
    noise = 0.5 * np.random.randn(len(spread_chips))
    received = spread_chips + jammer + noise

    # 逆拡散
    pn_full = np.tile(pn, N_bits)
    despread = received * pn_full

    # ビット判定
    detected = []
    for i in range(N_bits):
        start = i * N_chips
        end = (i+1) * N_chips
        integral = np.sum(despread[start:end])
        detected.append(1 if integral > 0 else -1)

    ber = np.mean(np.array(detected) != data_dsss)
    ber_dsss.append(ber)

# === FH-SS ===
M_slots = 31  # 周波数スロット数（処理利得と揃える）
data_fhss = 2 * np.random.randint(0, 2, N_bits) - 1
hop_pattern_fhss = np.random.randint(0, M_slots, N_bits)

# ジャマーは1つのスロットを妨害
jammer_slot = 5

ber_fhss = []
for js_db in jammer_power_range:
    js_lin = 10**(js_db/10)
    noise_base = 0.5

    detected = []
    for i in range(N_bits):
        # 受信信号 = data + noise + (ジャミング if same slot)
        rx = data_fhss[i] + noise_base * np.random.randn()
        if hop_pattern_fhss[i] == jammer_slot:
            rx += np.sqrt(2 * js_lin) * (2*np.random.randint(0,2)-1)
        detected.append(1 if rx > 0 else -1)

    ber = np.mean(np.array(detected) != data_fhss)
    ber_fhss.append(ber)

# パーシャルバンドジャミング（複数スロット妨害）
ber_fhss_partial = []
n_jammed_slots = 5  # 5/31スロットを妨害
jammed_slots = set(range(n_jammed_slots))

for js_db in jammer_power_range:
    js_per_slot = 10**(js_db/10) / n_jammed_slots  # 電力分散
    noise_base = 0.5

    detected = []
    for i in range(N_bits):
        rx = data_fhss[i] + noise_base * np.random.randn()
        if hop_pattern_fhss[i] in jammed_slots:
            rx += np.sqrt(2 * js_per_slot) * (2*np.random.randint(0,2)-1)
        detected.append(1 if rx > 0 else -1)

    ber = np.mean(np.array(detected) != data_fhss)
    ber_fhss_partial.append(ber)

# 結果プロット
fig, ax = plt.subplots(figsize=(10, 6))
ax.semilogy(jammer_power_range, np.array(ber_dsss) + 1e-4, 'bo-', label='DS-SS (narrowband jammer)')
ax.semilogy(jammer_power_range, np.array(ber_fhss) + 1e-4, 'r^-', label='FH-SS (single-slot jammer)')
ax.semilogy(jammer_power_range, np.array(ber_fhss_partial) + 1e-4, 'gs-', label=f'FH-SS (partial-band: {n_jammed_slots}/{M_slots} slots)')
ax.set_xlabel('Jammer-to-Signal Ratio [dB]')
ax.set_ylabel('Bit Error Rate')
ax.set_title('DS-SS vs FH-SS: Anti-jamming Performance')
ax.legend()
ax.grid(True, which='both', alpha=0.3)
ax.set_ylim([1e-4, 1])

plt.tight_layout()
plt.savefig('dsss_vs_fhss_jamming.png', dpi=150, bbox_inches='tight')
plt.show()