M序列

  M 序列是一种伪随机序列，具有很好的伪噪声特性，常用于信道噪声测试和保密通信。不过 M 序列还有一个用途，也就是本文所介绍的——通过 M 序列测量频率响应。在讨论这个问题之前，我们先介绍 M 序列的特征与生成方法。

  M 序列通过 N 阶线性反馈移存器生成，其周期为 $2^N-1$ 。

如上图所示，反馈到最左端的信号为
$$a_n=\sum _{i=1}^n{c}_i{a}_{n-i}(mod2)$$
经过 $k$ 次移位后的输出为
$$a_k=\sum _{i=1}^n{c}_i{a}_{k-i}(mod2)$$
  用特征多项式表示反馈的连接状态
$$f(x)=c_n{x}_n+c_{n-1}{x}^{n-1}+\cdots +c_{1}x+c_0=\sum _{i=0}^n{c}_i{x}^i$$
由于连接状态中，$c_0、c_n$ 均为 1，因此特征多项式必然包含 $x^n+1$ ，特征多项式的次数 n 即是移存器的级数。

  输出序列用多项式表示（母函数）
$$G(x)=a_0+a_{1}x+a_2{x}^2+\cdots =\sum _{i=0}^{\infty }{a}_i{x}^i$$
本原多项式：若一个多项式满足：1. $f(x)$ 是既约的；2. $f(x)$ 可以整除 $2^m+1$ ，$m=2^n-1$ ；3. $f(x)$ 除不尽 $2^q+1$ ，$q<m$ ；则称该多项式为本原多项式。

  当特征多项式为本原多项式时，反馈移存器的输出序列为 M 序列。

matlab 实现

% M 序列生成器
% @mg: M序列
% @f : 本原多项式
% a(k)=mod(sum_1^n(c(i)*a(k-i)),2)
function mg=m_generator(f)
    f=dec2bin(f)-'0';
    n=length(f)-1;    %M序列阶数
    N=2^n-1;          %伪随机码的周期
    register=[zeros(1,n-1),1];
    newregister=zeros(1,n);
    mg=zeros(1,N);
    f_inv=f(end-1:-1:1);
    for i=1:N
        newregister(1)=mod(sum(f_inv.*register),2);
        newregister(2:n)=register(1:n-1);
        register=newregister;
        mg(i)=register(n);
    end
end

% 测试
clc,clear,close all

N=10;   % N级反馈移存器，需要N-1阶本原多项式
L=2^N-1;
F=primpoly(N,'all'); %显示该阶数下的全部本原多项式

mg=m_generator(F(1)); %使用一个本原多项式生成m序列

M序列测量滤波器频率响应

  通常测量某一系统的幅频响应，是通过扫频来完成的，然而需要时间较长，还需要实现 DDS，比较复杂。而考虑到 m 序列具有白谱，因此测量经过滤波器后的 m 序列的值，并做 DFT，即可获知滤波器的幅频特性。

  假设生成 10 阶 M 序列（1023 点），以 $256kHz$ 的比特速率发送该序列，经过滤波器后，以同样的采样率用 ADC 采样，则频率分辨率为 $256kHz/1023=\sim 250Hz$ 。

%---------m序列--------
clc,clear,close all

N=10;   % N级反馈移存器，需要N-1阶本原多项式
L=2^N-1;
F=primpoly(N,'all'); %显示该阶数下的全部本原多项式

mg=m_generator(F(1))'*2-1;

%% m序列幅频特性
fs=256e3;
freq=fs*(0:(L/2))/L;

mg_fft=fft(mg)/L*2;
plot(freq/1e3,abs(mg_fft(1:floor(L/2)+1)))  %幅频显示为白谱
xlabel('f/kHz')

%% 滤波
mg2=smooth(mg);

mg2_fft=fft(mg2)/L*2;
figure
plot(freq/1e3,abs(mg2_fft(1:floor(L/2)+1))) %滤波器幅频特性曲线
xlabel('f/kHz')

测试中使用了平滑滤波器，测得的滤波器频率响应如下