python3で『PSK変調』と『ダイレクトコンバージョン』

 以下『Phase Shift Keying』のサンプルコードです。

誰かの役に立てば幸いなり。

＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃＃

import numpy as np

from scipy import signal as sg

from scipy.fftpack import fft,ifft,fftshift,fftfreq

import matplotlib.pyplot as plt

fb = 440 #baseband

fc = 10000 #carrier

fs = 44100 #sampling freq

A = 1 #amp

t = 3 #time

alpha = np.pi/4 

beta = np.pi/6

dat = np.arange(0,t,1/fs)

pt = 2*np.pi*dat

#sin(alpha)cos(beta)

sig = A * np.sin(fb * pt + alpha)

carrier = A * np.cos(fc * pt + beta)

# modulated signal

wave = A * np.cos((fc*pt) + beta + sig)

#plt.plot(sig[0:100])

#plt.plot(carrier[0:100])

#plt.plot(wave[0:100])

#plt.show()

#AWGN noise add 

#Add AWGN noise to the transmitted signal

nMean = 0 #noise mean

nSigma = 0.1 #noise sigma

n = np.random.normal(nMean, nSigma, len(dat))

noisy = wave + n  #noisy received signal

#plt.plot(wave[0:100])

#plt.plot(noisy[0:100])

#plt.show()

#analytic signal

asig = sg.hilbert(noisy)

# instaneous phase

phase = np.unwrap(np.angle(asig))

#demodulation

#

offset = 2 * np.pi * fc * dat + beta

want = phase - offset

#noise reduction

#low pass filter

# num must be odd number. num up more smooth

num = 7

nyq = fs/2

cutoff = (fc/2) / nyq

lpf = sg.firwin(num,cutoff)

signal = sg.lfilter(lpf,1,want)

plt.plot(signal[:300])

plt.plot(sig[:300])

plt.show()