Chirp信号公式与对离散生成算法之间的差异matlabTsinghuaJoking-

讨论产生线性频率变化的公式和它的离散公式之间的差异，并提出Chirp信号的修改方案。

Chirp信号的公式

对于固定频率
$f_1$

f1​的信号，它的表达式为：
$rleft( t right) = cos left( {2pi cdot f_1 cdot t} right)$

r(t)=cos(2π⋅f1​⋅t)

其中正弦信号中的相位是频率对时间的积分：
$theta left( t right) = 2pi cdot int_0^t {f_1 cdot dtau } = 2pi cdot f_1 cdot t$

θ(t)=2π⋅∫0t​f1​⋅dτ=2π⋅f1​⋅t

Chirp信号的频率是变化的，一种最简单的变化就是线性变化，从开始
$t = t_0$

t=t0​的频率
$f_{start}$

fstart​，一直线性变化到
$t = t_1$

t=t1​时的频率
$f_{end}$

fend​。那么任何时刻的频率值为：
$fleft( t right) = {{f_{end} – f_{start} } over {t_1 – t_0 }} cdot t + f_{start}$

f(t)=t1​−t0​fend​−fstart​​⋅t+fstart​

对应的正弦信号相角为：
$theta left( t right) = 2pi cdot int_{t_0 }^{t_1 } {fleft( tau right)dtau } = 2pi int_{t_0 }^{t_1 } {{{f_{end} – f_{start} } over {t_1 – t_0 }} cdot tau + f_{start} dtau }$

θ(t)=2π⋅∫t0​t1​​f(τ)dτ=2π∫t0​t1​​t1​−t0​fend​−fstart​​⋅τ+fstart​dτ
$= 2pi left[ {{{f_{end} – f_{start} } over {t_1 – t_0 }} cdot {1 over 2}left( {t_{}^2 – t_0^2 } right) + f_{start} left( {t_{} – t_0 } right)} right]$

=2π[t1​−t0​fend​−fstart​​⋅21​(t2​−t02​)+fstart​(t​−t0​)]

如果假设
$t_0 = 0$

t0​=0，那么上面的公式就是：
$theta left( t right) = 2pi left[ {{{f_{end} – f_{start} } over {t_1 }} cdot {1 over 2}t^2 + f_{start} cdot t} right]$

θ(t)=2π[t1​fend​−fstart​​⋅21​t2+fstart​⋅t]

所以，Chirp信号的公式为：

$rleft( t right) = cos left[ {2pi cdot left( {{{f_{end} – f_{start} } over {t_1 }} cdot {1 over 2}t^2 + f_{start} cdot t} right)} right]$

r(t)=cos[2π⋅(t1​fend​−fstart​​⋅21​t2+fstart​⋅t)]

下面是在智能车比赛中相应的参数：

起始时间：
$t_0 = 0s$

t0​=0s
结束时间：
$t_1 = 0.2048s$

t1​=0.2048s
起始频率：
$f_{start} = 250Hz$

fstart​=250Hz
结束频率：
$f_{end} = 2000Hz$

fend​=2000Hz

利用时间间隔
$t_s = 1/f_s = 100mu s$

ts​=1/fs​=100μs对该信号进行离散化，得到的波形为：

def chirpf(t): return cos(2*pi*((fend-fstart)/t1*0.5*t*t+fstart*t)) 

使用MATLAB产生的Chirp信号波形：

tend=0.2048 fs = 1e4 fstart=250 fend=2000 tdim = 0:1/fs:(tend-1/fs) chirpdim=chirp(tdim,fstart,tend,fend) 

import scipy.signal data = scipy.signal.chirp(tdim, fstart,t1, fend) 

通过数值对比，可以看到上面三种方式得到的Chirp信号的误差为0。下图显示了公式所计算得到的Chirp信号与MATLAB所产生的Chirp信号的差值。其中微小的误差可能来源于MATLAB结果经过clipboard剪切板之后，数值显示精度降低所引起的。

另外一部分则可能是对于数据的的最后一点的频率微小差异的定义。

Chirp信号离散生成算法

下面是智能车竞赛信号板单片机通过DAC输出Chirp信号对输出信号buffer进行初始化的程序。

//------------------------------------------------------------------------------ void InitDACBuffer(float fStartF, float fEndF) { float fAngle = 0; float fFrequency; float fDeltaT = 1.0 / DAC_OUTPUT_FREQUENCY; int i; for(i = 0; i < DAC_BUFFER; i ++) {                 g_nDACBuffer[i] = (unsigned short)((sin(fAngle * 2 * 3.1415926) + 1.0) / 2 * 0x4ff)+0x100;         fFrequency = (fEndF - fStartF) * (i + 1) / DAC_BUFFER + fStartF;         fAngle += fFrequency * fDeltaT; } } 

在程序初始化期间调用如下命令：

 InitDACBuffer(250, 2000); 

根据上面产生的程序，下面绘制出对应的波形：

DAC_BUFFER = 2048 dacbuf = [] angle = 0 deltat = 1/1e4 for i in range(0, DAC_BUFFER):     data = cos(angle*2*pi)     dacbuf.append(data)     frequency = (fend - fstart) * (i + 1) / DAC_BUFFER + fstart     angle = angle + frequency * deltat  plt.plot(tdim, dacbuf) plt.xlabel('Time(s)') plt.ylabel('Amplitude') plt.grid(True) plt.show() 

对比由离散程序所产生的数据与前面通过公式计算出的波形，可以看到它们之间存在着差异。下面将公式所计算出的数值与前面程序所产生信号相加，所得到的误差信号：

使用程序产生Chirp信号的时候，信号的相位是通过离散积分而得到的，其中所出现的积分误差随着时间增加而叠加，最终使得所产生的信号的相位与实际的相位越来越大，信号的差别，也就增加了。

修改方案

将最终信号Chirp生成的公式修改为：

$t_n = n cdot t_s = {n over {f_s }} = n cdot 10^{ – 4}$

tn​=n⋅ts​=fs​n​=n⋅10−4

$xleft[ n right] = sin left[ {2pi cdot left( {{{f_{end} – f_{start} } over {t_1 }} cdot {1 over 2}t_n^2 + f_{start} cdot t_n } right)} right]$

x[n]=sin[2π⋅(t1​fend​−fstart​​⋅21​tn2​+fstart​⋅tn​)]

使用sin信号的目的，就是防止最初的信号开始的跳变。根据上面公式所产生的波形为：

生成12bit DAC转换对应的数据：

for i in range(DAC_BUFFER):     tn = i / fs     angle = (fend-fstart)/t1 * tn * tn / 2 + fstart * tn     datasin = sin(2*pi*angle)      dataint16 = int((datasin + 1.0) / 2 * 0x4ff + 0x100)     dacbuf.append(dataint16) 

信号波形为：

生成相关的整型信号，便于MCU输出模拟信号：

#------------------------------------------------------------ fstart = 250 fend = 2000 t0 = 0 t1 = 0.2048 ts = 100e-6  tdim = linspace(t0, t1, int((t1-t0)/ts), endpoint=False) #------------------------------------------------------------ DAC_BUFFER = 2048 dacbuf = [] fs = 1e4 for i in range(DAC_BUFFER):     tn = i / fs     angle = (fend-fstart)/t1 * tn * tn / 2 + fstart * tn     datasin = sin(2*pi*angle)      dataint16 = int((datasin + 1.0) / 2 * 0x3f)     dacbuf.append(dataint16)  tspsave('chirp6bit', dacbuf=dacbuf) 

下面是将NPZ数据转换成C51的变量的程序。

dacbuf0 = tspload('chirp8bit', 'dacbuf') dacbuf1 = tspload('chirp7bit', 'dacbuf') dacbuf2 = tspload('chirp6bit', 'dacbuf')  pastestr = 'unsigned char const code g_ucChirpData[] = {rn' linenumber = 16 lines = int(len(dacbuf2) / linenumber) for i in range(lines):     linedata = 127 - dacbuf2[i * linenumber : (i + 1) * linenumber]      datastr = str(linedata).strip('[]').split(' ')     datastr1 = [s for s in datastr if len(s) > 0]     linestr = ','.join(datastr1).lstrip(',') + ',rn'        pastestr = pastestr + '    ' + linestr  pastestr = pastestr + '};rn'  clipboard.copy(pastestr) printf(pastestr)