近期在信号与系统课程中讲完了“信号的采样与恢复”的内容。通常情况下对于信号的采样都是沿着时间轴对信号的幅值进行采样，获得信号的离散时间点上的数据。
如果将信号的波形绘制在直接坐标系中，那么该曲线就是分布在二维空间上的曲线。曲线上的点可以沿着时间轴进行排列，当然也可以按照幅值的大小进行排列。也就是按照取值间隔，将信号通过该间隔的时间进行保存，这就是对信号的时间采样。
沿着时间轴对信号的幅值进行采样，Nyquist-Shannon定理告诉我们如何进行采样和如何进行信号恢复。但是，如何对信号时间进行采样，如何恢复是一个经典的未解决问题。
Logan定理[Logan,Jr.，1977]对一种特殊信号给出了采样与恢复的描述：如果一个信号的频谱具有倍频的性质，即信号频谱分布在一个频率范围内，最高频率是最低频率的两倍。那么这个信号可以通过它的过零点的时间值进行恢复。恢复的信号与原始信号仅仅相差一个比例因子。
信号采样与恢复
1. 幅度采样
下面是一段普通语音信号的复制采样波形。通过简单的DA转换和低通滤波便可以恢复出原始的声音波形。
这段语音信号的内容为：
2. 时间采样
如果仅仅保留该信号的过零时间点的信息，它的幅值全部去掉，所形成的波形大体上如下图所示。这很像将原来的语音信号通过一个过零点比较器，输出的信号反映了信号的极性。
▲ 声音信号过零点采样
当然，这个信号中包含有原来信号的部分信息。但它并不是原来语音信号的完美的回复。
尽管如此，播放这个信号，还是可以听到原来语音的信息。虽然有很大的失真和噪声。这说明原来的信号信息还是部分保留在这些过零点中。
经过比较器之后的语音信号为：
如果语音信号满足Logan定理的要求，那么理论上是可以恢复出原来的信号的。但如何来恢复？
如果信号本身是一个周期信号，也就是信号的频谱是离散的频谱。相对回复信号的算法比较简单。Sam Roweis等人在论文“Signal Reconstruction from Zero-Crossings"中给出了通过求解数据矩阵零空间向量的方法，来通过信号的过零时间点来重构信号的方法。
重建算法
1. 基本原理
已知到信号具有倍频窄带频谱，它的频率范围分布在 范围内。
已知信号的周期，，以及个信号过零点：
重构信号的计算步骤如下：
(1)计算相关参数：以及 ；
(2)构造矩阵：；
(3)寻找数据矩阵零空间向量：构造数据矩阵，矩阵大小是。对该矩阵进行奇异值分解，得到。其中是空间上的正交矩阵，是空间上的正交矩阵。是奇异值向量，长度为。
零空间向量是奇异值向量中最小(理论上应该为0，但由于计算误差的存在，它可能是一个很小的数)对应的中的向量，由于SVD算法往往把结果按照绝对值从大到小排列，所以的最后一个向量就对应着数据矩阵的零空间向量。
(4)构造信号函数：将数据空间矩阵中零空间向量前个数值当做，后个数值当做，重建信号公式为：
这种回复信号的过程，实际上就是根据信号的过零点来求解上面的函数中的参数。具体的理论分析在这里就不再展开了。
2. 测试函数
(1)实验信号的数学表达式：
选择一个频率分布在10Hz到20Hz之间的一个信号进行实验，随机指定对应的cos,sin信号的系数，如下：
这是一个周期为1的倍频信号。
(2)信号的产生Python程序：
使用下面python程序，可以产生该信号的数据。也可以通过该函数完搜索信号的过零点。
def sfunc1(x): pi2 = 2 * pi retdata = cos(pi2*11*x) + sin(pi2*11*x)/2 + \ cos(pi2*12*x)/33 + sin(pi2*12*x)/4 + \ cos(pi2*13*x) + sin(pi2*13*x)/8 + \ cos(pi2*14*x) + sin(pi2*14*x)/7 + \ cos(pi2*15*x)/22 + sin(pi2*15*x)/3 +\ cos(pi2*16*x) + sin(pi2*16*x)/12 +\ cos(pi2*17*x) + sin(pi2*17*x)/40 +\ cos(pi2*18*x) + sin(pi2*18*x)/2 +\ cos(pi2*19*x)/3 + sin(pi2*19*x)/2 return retdata
(3)实验信号的波形：
下面绘制出0~2秒两个周期内的波形。
绘制该信号的过零点饱和信号，它仅仅保留了该信号的过零点的时间和相位信息。计算公式为：
(4)信号的过零点：
通过数值计算，来获得信号的过零点。下面重新绘制出信号一个周期内的波形，没有添加任何噪声。
通过对区间(0,1)采集10^6^个数值，然后通过寻找过零点，获得二十八个信号的过零点的值。
搜寻函数值过零点的python程序如下crosszero(t,val)。其中是函数的自变量，是函数值的采样。函数返回是对应函数过零点时的的数值。
def crosszero(t, val): valsign = sign(val) valsignchange = [int(x!=y) for x,y in zip(valsign[0:-1],valsign[1:])] tvalue = [(x,y) for x,y in zip(t[0:-1], valsignchange)] zerot = filter(lambda t: t[1]!= 0, tvalue) return [zt[0] for zt in zerot]
通过scipy.optimize.root来寻找信号的根，用于确定信号的过零点。利用上面搜索的结果作为初始值。
sol = scipy.optimize.root(tssub.sfunc1, czt, method='lm')
tssub.sfunc1：定义的信号函数；czt：是前面通过数值过零点搜索获得的28个根的数值；
如下是sol['x']中的数值，包含了最终优化后的数值，对比前面通过搜索获得数值，可以看到基本上在10^-6^的内存在一定的误差。
将通过数值计算所得到的28个函数的根绘制在信号波形上，看到他们的分布。
3. 重建结果
根据前面所叙述的方法，使用28个过零点信息重构出的信号波形如下图所示。重构的信号与原始的信号之间波形基本一致，只是相差了一个比例因子。
前面实验中的程序和数据可以在CSDN博文中看到。
http://zhuoqing.blog.csdn/net
