numpy  ·  22日 7月 2020

    パルスで直交復調(5) -numpyでやってみる-

    波形生成

    fsk_study_gen_ver2.py

    1. # 変数クリア
    2. from IPython import get_ipython
    3. get_ipython().magic('reset -sf')
    4.  
    5. import matplotlib.pyplot as plt
    6. import numpy as np
    7. import scipy.signal as signal
    8.  
    9. # シンボル列を引き延ばす
    10. def extend_symbol_to_simulation_sampling_rate(symbols,t_symbol,t_simulation_sampling):
    11.     # シンボル数
    12.     size_of_symbols=symbols.size
    13.     # 最後のシンボルが終わる時間
    14.     end_time=size_of_symbols*t_symbol
    15.     # 計算タイミング
    16.     t=np.arange(0,end_time,t_simulation_sampling)
    17.     # 各シンボルの開始（終了）時間
    18.     t_periods_of_symbols=np.append(0,np.arange(1,size_of_symbols,1)*t_symbol)
    19.     # 結果データ領域確保
    20.     ex_symbols=np.zeros_like(t)
    21.     # 先頭データは入れておく
    22.     ex_symbols[0]=symbols[0]
    23.     # 各シンボルについて、その期間のインデックスを取得し、データを入れ込む
    24.     for i in range(0,size_of_symbols):
    25.         ex_symbols[np.where(t>t_periods_of_symbols[i])]=symbols[i]
    26.     return ex_symbols,end_time
    27.          
    28. # よく使う変数
    29. pi=np.pi
    30. deg2rad=pi/180.0
    31. twopi=2*pi
    32.  
    33. # 設定
    34. f_carrier=433.92e6              # 搬送波周波数 in Hz
    35. f_sym=1e6                       # シンボルレート in Hz
    36. m=0.5                           # 変調指数
    37. f_samp_simulation=f_carrier*10  # 計算のサンプリングレート in Hz
    38. phase_initial=0*deg2rad         # 搬送波の初期位相(何でもいい) in rad
    39. f_deviation=(m*f_sym)/2         # 周波数偏差 in Hz
    40.  
    41. symbols_in_base=signal.max_len_seq(5)[0] # シンボル
    42. #symbols_in_base[0::2]=1
    43. #symbols_in_base[1::2]=0
    44.  
    45. # 後で使う変数
    46. t_sym=1/f_sym
    47. t_samp_simulation=1/f_samp_simulation
    48. omega_carrier=twopi*f_carrier
    49. omega_deviation=twopi*f_deviation
    50.  
    51. # シンボルを計算用に拡張
    52. symbols,t_end=extend_symbol_to_simulation_sampling_rate(symbols_in_base,t_sym,t_samp_simulation)
    53. symbols=symbols*2-1 # -1 or +1 にする
    54.  
    55. # 計算するタイミング
    56. t=np.arange(0,t_end,t_samp_simulation)
    57.  
    58. # 全計算タイミングにおける瞬間的な角速度
    59. omega_in_moment=omega_carrier+omega_deviation*symbols
    60.  
    61. # 全計算タイミングにおける瞬間的な位相変化量
    62. phase_change_in_moment=omega_in_moment*t_samp_simulation
    63.  
    64. # 全計算タイミングでの位相
    65. phase_t=np.cumsum(phase_change_in_moment)+phase_initial
    66.  
    67. # 波形生成
    68. amp_t=np.exp(1j*phase_t)
    69.  
    70. np.save('output/fsk_o_t',t)
    71. np.save('output/fsk_o_amp_t',amp_t)
    72. np.save('output/fsk_o_f_sym',f_sym)
    73. np.save('output/fsk_o_f_carrier',f_carrier)
    74.  
    75. # プロット
    76. fig=plt.figure()
    77. fig.add_subplot(2,1,1)
    78. plt.plot(t,symbols)
    79. fig.add_subplot(2,1,2)
    80. plt.plot(t,amp_t)

     

    普通に正弦波の局発でやってみる

    fsk_study_qdemod_ver2.py

    1. # 変数クリア
    2. from IPython import get_ipython
    3. get_ipython().magic('reset -sf')
    4.  
    5. import matplotlib.pyplot as plt
    6. import numpy as np
    7. from scipy import signal
    8.  
    9. # 波形を読み込む
    10. f_carrier=np.load("output/fsk_o_f_carrier.npy")
    11. f_sym=np.load("output/fsk_o_f_sym.npy")
    12. t=np.load('output/fsk_o_t.npy')
    13. amp_complex_t=np.load('output/fsk_o_amp_t.npy')
    14. t_sym=1/f_sym
    15. amp_t=amp_complex_t.real
    16. td0=t[:-1:]
    17. td1=t[1::]
    18. t_samp_simulation=np.average(td1-td0)
    19. f_samp_simulation=1/t_samp_simulation
    20.  
    21. # よく使う変数
    22. pi=np.pi
    23. deg2rad=pi/180.0
    24. twopi=2*pi
    25.  
    26. # 設定
    27. f_Lo=f_carrier            # 搬送波周波数 in Hz
    28. phase_Lo=0*deg2rad    # 搬送波の初期位相(何でもいい) in rad
    29. m=0.5 # 波形生成時の変調指数
    30. f_deviation=(m*f_sym)/2
    31.  
    32. # 後で使う変数
    33. t_samp_simulation=1/f_samp_simulation
    34. omega_Lo=twopi*f_Lo
    35.  
    36. # 局発信号の生成
    37. amp_LoI_t=np.cos(omega_Lo*t+phase_Lo)
    38. amp_LoQ_t=np.sin(omega_Lo*t+phase_Lo)
    39.  
    40. # ミキサー
    41. amp_I_t_before_LPF=amp_t*amp_LoI_t
    42. amp_Q_t_before_LPF=amp_t*amp_LoQ_t
    43.  
    44. # フィルター係数生成
    45. freq_cutoff=f_sym*2
    46. w_cutoff=freq_cutoff/(f_samp_simulation/2)
    47. b,a=signal.butter(1,w_cutoff,'lowpass')
    48.  
    49. # フィルター適用
    50. amp_I_t=signal.lfilter(b,a,amp_I_t_before_LPF)
    51. amp_Q_t=signal.lfilter(b,a,amp_Q_t_before_LPF)
    52. amp_demod_t=amp_I_t+1j*amp_Q_t
    53.  
    54. # ダウンサンプリング
    55. f_sampling_target=f_sym*10
    56. downsampling_rate=int(f_samp_simulation/f_sampling_target)
    57. t_downsample=t[::downsampling_rate]
    58. amp_demod_t_downsample=amp_demod_t[::downsampling_rate]
    59. f_samp_simulation_downsample=f_samp_simulation/downsampling_rate
    60. t_samp_simulation_downsample=1/f_samp_simulation_downsample
    61.  
    62. # 復調
    63. amp_demod_t0_downsample=amp_demod_t_downsample[:-1:]
    64. amp_demod_t1_downsample=amp_demod_t_downsample[1::]
    65. phase_shift_t_downsample=-1*np.angle(amp_demod_t1_downsample/amp_demod_t0_downsample)
    66. phase_shift_t_downsample=np.append(0,phase_shift_t_downsample)
    67.  
    68. # 規格化
    69. phase_shift_t_downsample=phase_shift_t_downsample/(f_sym/f_samp_simulation_downsample*pi*m)
    70.  
    71. # プロット
    72. fig=plt.figure()
    73. fig.add_subplot(2,1,1)
    74. plt.plot(t_downsample,amp_demod_t_downsample.real)
    75. plt.plot(t_downsample,amp_demod_t_downsample.imag)
    76. fig.add_subplot(2,1,2)
    77. plt.plot(t_downsample,phase_shift_t_downsample)
    78. plt.ylim(-2,2)
    79.  
    80. fig=plt.figure()
    81. ax=fig.add_subplot(1,1,1)
    82. plt.scatter(amp_demod_t_downsample.real,amp_demod_t_downsample.imag)
    83. plt.xlim(-0.5,0.5)
    84. plt.ylim(-0.5,0.5)
    85. ax.set_aspect('equal', adjustable='box')

     

     

    矩形波でやってみる。

    fsk_study_qdemod_lo_sqwave_ver2.py

    1. # 変数クリア
    2. from IPython import get_ipython
    3. get_ipython().magic('reset -sf')
    4.  
    5. import matplotlib.pyplot as plt
    6. import numpy as np
    7. from scipy import signal
    8.  
    9. # 波形を読み込む
    10. f_carrier=np.load("output/fsk_o_f_carrier.npy")
    11. f_sym=np.load("output/fsk_o_f_sym.npy")
    12. t=np.load('output/fsk_o_t.npy')
    13. amp_complex_t=np.load('output/fsk_o_amp_t.npy')
    14. t_sym=1/f_sym
    15. amp_t=amp_complex_t.real
    16. td0=t[:-1:]
    17. td1=t[1::]
    18. t_samp_simulation=np.average(td1-td0)
    19. f_samp_simulation=1/t_samp_simulation
    20.  
    21. # よく使う変数
    22. pi=np.pi
    23. deg2rad=pi/180.0
    24. twopi=2*pi
    25.  
    26. # 設定
    27. f_Lo=f_carrier          # 搬送波周波数 in Hz
    28. phase_Lo=0*deg2rad      # 搬送波の初期位相(何でもいい) in rad
    29. m=0.5                   # 波形生成時の変調指数
    30. f_deviation=(m*f_sym)/2
    31.  
    32. # 後で使う変数
    33. t_samp_simulation=1/f_samp_simulation
    34. omega_Lo=twopi*f_Lo
    35.  
    36. # 局発信号の生成
    37. amp_LoI_t=np.cos(omega_Lo*t+phase_Lo)
    38. amp_LoQ_t=np.sin(omega_Lo*t+phase_Lo)
    39. amp_LoI_t[amp_LoI_t<0]=-1/4*pi
    40. amp_LoI_t[amp_LoI_t>=0]=1/4*pi
    41. amp_LoQ_t[amp_LoQ_t<0]=-1/4*pi
    42. amp_LoQ_t[amp_LoQ_t>=0]=1/4*pi
    43.  
    44. # ミキサー
    45. amp_I_t_before_LPF=amp_t*amp_LoI_t
    46. amp_Q_t_before_LPF=amp_t*amp_LoQ_t
    47.  
    48. # フィルター係数生成
    49. freq_cutoff=f_sym*2
    50. w_cutoff=freq_cutoff/(f_samp_simulation/2)
    51. b,a=signal.butter(1,w_cutoff,'lowpass')
    52.  
    53. # フィルター適用
    54. amp_I_t=signal.lfilter(b,a,amp_I_t_before_LPF)
    55. amp_Q_t=signal.lfilter(b,a,amp_Q_t_before_LPF)
    56. amp_demod_t=amp_I_t+1j*amp_Q_t
    57.  
    58. # ダウンサンプリング
    59. f_sampling_target=f_sym*10
    60. downsampling_rate=int(f_samp_simulation/f_sampling_target)
    61. t_downsample=t[::downsampling_rate]
    62. amp_demod_t_downsample=amp_demod_t[::downsampling_rate]
    63. f_samp_simulation_downsample=f_samp_simulation/downsampling_rate
    64. t_samp_simulation_downsample=1/f_samp_simulation_downsample
    65.  
    66. # 復調
    67. amp_demod_t0_downsample=amp_demod_t_downsample[:-1:]
    68. amp_demod_t1_downsample=amp_demod_t_downsample[1::]
    69. phase_shift_t_downsample=-1*np.angle(amp_demod_t1_downsample/amp_demod_t0_downsample)
    70. phase_shift_t_downsample=np.append(0,phase_shift_t_downsample)
    71.  
    72. # 規格化
    73. phase_shift_t_downsample=phase_shift_t_downsample/(f_sym/f_samp_simulation_downsample*pi*m)
    74.  
    75. # プロット
    76. fig=plt.figure()
    77. fig.add_subplot(2,1,1)
    78. plt.plot(t_downsample,amp_demod_t_downsample.real)
    79. plt.plot(t_downsample,amp_demod_t_downsample.imag)
    80. fig.add_subplot(2,1,2)
    81. plt.plot(t_downsample,phase_shift_t_downsample)
    82. plt.ylim(-2,2)
    83.  
    84. fig=plt.figure()
    85. ax=fig.add_subplot(1,1,1)
    86. plt.scatter(amp_demod_t_downsample.real,amp_demod_t_downsample.imag)
    87. plt.xlim(-0.5,0.5)
    88. plt.ylim(-0.5,0.5)
    89. ax.set_aspect('equal', adjustable='box')

     

    良いっちゃ～良いんだけど、、、やっぱり正弦波よりは美的に劣るか、、、

    ていうか、今までやったのパルスじゃなくて矩形波じゃん；

    カテゴリ： 非現実