2023全国大学生电子设计竞赛H题全解 [原创www.cnblogs.com/helesheng]

2023年又是全国大学生电子设计竞赛年，一如既往的指导学生死磕H题。8月2日看到公布的赛题，我自己还沾沾自喜，觉得今年学生用嵌入式系统和数字信号处理知识就可以完成这题，赛前都辅导过，应该成绩不差。哪想到结果大跌眼镜，不但成绩还不如往年，学生的解题思路更是各式各样……

直到这几天寒假才有空，我按照自己的思路， 在STM32上用数字信号处理的方法将这道题全部做了一遍，效果不错。现分享出来以飨读者和粉丝，博大家一笑尔！也供后续的学生和参赛者参考。（注：以下讲解未按照电赛组委会要求的报告标准格式撰写，而是使用了随笔形式，重点讲解结题思路。未来阅读此文的参赛者不可以此作为《设计报告》范本）

作为电子设计竞赛的省级评委以及一个从业多年教育工作者，在竞赛评奖中，我一方面看到了很多兄弟学校参赛队对H题的不同解法，获益良多。另一方面，我也发现了很多高校在电子电工、信号处理前、嵌入式系统等领域的教学短板。因此我会在接下来的解题过程中，在每个关键步骤之后总结一下相关教学中存在的问题，希望对我省的电子电工教学水平的提高有所帮助。

先“叠个甲”：本文观点难免有失偏颇，有不妥之处还请各位教育专家大佬不要介怀，并批评指正。以下原创内容欢迎网友转载，但请注明出处：https://www.cnblogs.com/helesheng 

 

1.       2023全国大学生电子设计竞赛H题原题

 

下面是原题：

2.       审题和选择技术路径

 本题的核心功能是题图1中的分离电路实现的信号 “分离”，光从名称理解，很容易让人想到对C信号进行滤波——把两个信号分开。但仔细思考后会发现，这样做的问题是滤波器的过渡带极窄，极端情况下95KHz和100KHz之间只相差5%。很难将两个信号分离干净：A'中含有B的频率；B'中含有A的频率。从而造成A'和B'两个信号的失真。提高信号分离度的有效方法只有提高滤波器的阶数：对模拟滤波器而言，由于电容的精度很难做到20%以上，高阶模拟滤波器很难实现准确的截止频率和很窄的过渡带；对数字滤波器而言，高阶的IIR滤波器很难稳定，高阶的FIR滤波器会使计算量大增，无法满足计算的实时性。（注：此处必须具备一定的信号处理实践知识，既要能够判断5%的过渡带很难实现，又要明白IIR的缺点在于稳定性，FIR滤波器的短板在于计算量是IIR的许多倍）

另一种实现思路是“重建”A'和B'两个信号，即让图1中的分离电路重新生成C信号中含有的A、B两个信号。有很多参赛者怀疑重建的信号是否满足题目要求。对此组委会官方在赛题公布后公开答疑，明确了重建信号是满足要求的。

“重建”信号方案的难度在于：由于信号发生器产生的A、B和重新产生的A'、B'信号都不可能是频率绝对准确的，其精度都在0.01%左右，也就是万分之一。对于100KHz的信号，最糟糕的情况下，两者之间每秒将相差±10个左右的波形。如果用示波器同时显示两个信号会出现如下现象：用A信号作为触发源，则A'信号每秒在屏幕上将向左或向右“滚动”最多20个波形，也就是无法做到原题第三部分“说明”第（5）条提到的用示波器能够“稳定不漂移”的显示原信号和重建信号。

解决问题的办法是：不断地实时测量A与A'、B与B'之间的相位差，并根据这个相位差不断调整重建的A'和B'的相位。由此可以产生多种不同的技术路径及相关的解决方案：

 1. 全数字方案：对C信号采样，在数字域计算A、B的频率，并用DAC重建同频率的两路信号A'、B'。随后不停地重复计算输入的A、B与A'、B'的相位差，并调整A'、B'的相位。

2. 将测量相位差的工作交给模拟相位测量芯片（如AD8302）来完成，MCU通过其输出来获取实时相位差。

3. 把信号重建工作交给DDS(Direct Digital Frequency Synthesis，如AD9833/AD9854等)芯片来完成，MCU通过接口实时调整DDS产生信号的相位。

4. 把信号重建工作交给压控振荡器（VCO，voltage-controlled oscillator)来实现，VCO控制电压由相位差测量电路来产生。（注：由于VCO很难产生100KHz以下的标准正弦波，本方案不太容易获得较好结果。）

非常高兴看到，本省的参赛者采用了上述除第4种以外的所有技术路径，有一个参赛队结合第2和第3种技术路径获得了相当不错的结果，获得了省级一等奖，而其他参赛队表现都不如人意。接下来我将用全数字方案来解决这个问题。具体使用的器件是我校教学中使用的STM32F1系列芯片，熟悉其他器件平台的读者可以参考我的算法，在你熟悉的平台上实现。算法对于平台的基本要求是：至少具有一路1MSPS以上的A/D转换，两路1MSPS以上的D/A转换能力，以及片上乘加运算（MAC）指令，处理能力达到50MIPS以上。这些软硬件能力，目前主流的国产或知名大厂的大部分32位和部分16位嵌入式MCU或DSP平台都具备。

本环节揭示的信号处理教学中的问题： 

a) 在测试中发现很多参赛队采用了滤波器，甚至模拟滤波器作为解决这道题的基本技术路径。这揭示出我们的《信号与系统》、《数字信号处理》和《模拟电路基础》等课程的教学停留在理论教学层面——学生对于滤波器的过渡带有认知，但没有确切的感性认识。只是因为上课时讲过这个内容，参赛同学就认为滤波器能够解决实际工程中的一切问题，缺乏解决“复杂工程问题”的能力。

b) 实验教学中，学生接触和使用示波器的机会不足。导致其不理解示波器的触发功能。对波形“滚动”的成因及对策没有认知。

3.       模拟加法电路

题图1中的加法器使用运算放大器实现，电路如下图所示。电路本身没有什么技术含量，就是教科书上的加法器电路。但由于MCU上集成的A/D只能接受正电压，这里用了隔直电容C1和C2去除输入信号的直流偏置，并用电阻R4和R5分压产生的Vcc/2作为加法器输出信号的直流偏置。这样不管输入信号A和B的直流偏置是0还是其他电压，只要幅度合适，都可以被MCU片上的A/D转换器采集。

图3 带直流偏置的模拟加法器

图3电路图左下角的是MCU处理后输出的DAC_CH1和DAC_CH2两个信号的LC低通滤波器，其作用是平滑D/A转换输出的阶梯。而整个模块使用了右上角的PMOD接口来和MCU板（图3未画出）相连。实物连接示意图如下图所示。

 图4 模拟加法器实物连接示意图

本环节揭示的电子电工教学中的问题：

a) 模拟电路教学过于重视分立元件教学，对于运算放大器这一最重要的教学对象关注度不足。其实现代模拟系统中，除电力电子专业之外，直接使用分离的二、三极管的场合已经非常有限，没有必要因循守旧，抱持几十年不变的教学内容。

b) 学生对于“模拟电路工作点”，这一重要概念理解不深刻，不会用。原因是模拟电路教材中多数使用双电源电路，工作点天生为GND，不需要调整。但对于现代低功耗、手持式电路而言，多数使用单电源供电，调整模拟电路的工作点就成为必须掌握的技能。对此，模拟电路教材必须有所改变。

4.  嵌入式程序的总体框架

本题的嵌入式程序要完成两件事：其一判断输入C波形中的含有的信号A、B的频率和类型（正弦或三角波）；其二根据实时输入的C信号，不断的调整重建信号的相位，以保持二者的相位差基本不变。根据题目中“说明”部分（3）给出的要求可知，判断输入信号中两个频率成分的频率和类型的工作应该在按下“启动键”后完成，且只需要执行一次，对应于下图左侧部分。而实时调整重建信号的相位的工作则应该在判断完成后不断地重复执行，以保证后续随时用示波器观测A和A’、B和B’信号时都可以做到“稳定显示不漂移”。其程序流程对应下图中右侧的循环部分。

 

图5 嵌入式程序整体流程图

5.  相位误差算法以及相位校正的嵌入式实现

1）相位误差算法理论分析

正如本文第2小节“审题”中分析的，本题最大的难点在于如何实时地计算算法重新生成单频的A'、B'信号和C中的A、B信号之间的相位差，并在随后的计算不断实时校正，将相位差保持在一个固定的值，从而实现在示波器上“稳定显示不漂移”。我通过对输入信号C进行高速A/D采样数字化后，再用信号处理算法计算C中的A、B和重建的A'和B'之间的相位差。

图6 理论上的C信号和合成它的A、B两个信号 

说到测量某个固定频率信号的相位，大部分人会首先想到FFT这个工具，原因是FFT算法是“数字信号处理”课程重点介绍的内容，且大部分MCU有FFT算法库可以直接调用。我最初也想象从FFT的相位谱还原A、B信号初始相位的方法：

在时域内截取指定长度（L）的C信号，并用FFT计算该信号的相位谱，然后在相位谱中读取A、B信号频率位置对应对应的相位值。

但不幸的是，我在MATLAB中仿真这一算法时，发现结果根本不对！分析频域信息时才发现了下图所示的现象。

图7 FFT后的功率谱（兰）和相位谱（红）

仔细想来，这种现象也不奇怪：当采样率为1MSPS，信号长度为1024时，频率分辨率（也就是频域中两个点的间隔）为fs/L = 1M/1024 = 976.56Hz。理论上说，100KHz的信号所对应的频谱位置将在第102.4个点上，而离散傅里叶变换只能得到整数点上的频域信息，那么100KHz的时域信号能量将泄漏到第102和103个频谱位置上。这个问题对于功率谱并不重要，因为功率是标量可以直接对周边几个点上的功率谱直接叠加即可；但对于相位谱而言再想从周边几个点的相位谱恢复出非整数点对应频率（如100KHz）的相位就不容易了……

其实，多年前我还是学生的时候，曾经仔细思考过这个问题，后来看了胡广书老师的经典教科书《数字信号处理——理论算法与实现》^[1]第三章第八小节《对正弦信号截断的原则》才对此释然的。胡老师的教科书是这样写的“抽样频率应为信号频率的整数倍，抽样点数应包含整周期”。看到胡老师的书，又涌起“身不能至，心向往之”的崇拜之情。

图8 《数字信号处理——理论算法与实现》第三章第八小节《对正弦信号截断的原则》

至此若仍然坚持FFT方案就需要把采样长度L时间内输入信号的周期数调整为K个。由于采样率fs，截取长度L，截取时长内周期数K和A、B信号频率f之间必须满足下式：

 要保持K为整数就要保证（1）式的右边也是整数。但不幸的是为执行FFT快速算法，L一般应该是2的整数次幂；f是题目中给出的10KHz、15KHz、20KHz……95KHz和100KHz等固定值。也就要求采样率fs是一个非常精确的小数。如：L为1024，f为100KHz，K为100时fs为1.024MSPS；L为1024，f为95KHz，K为100时fs为0.9728MSPS。这种连续非整数的采样率fs，对于MCU集成的逐次逼近SAR式 A/D转换器几乎是无法做到的。（注：以常见的STM32F1系列MCU为例，在56MHz主频下，其片上集成的SAR A/D转换器使用DMA+硬件连续触发模式，可获得最高采样率1MSPS。改变MCU主频、A/D工作频率或采样触发模式可获得不同的采样率，但在100KSPS以上的高采样率下，无法再获准确的获得连续变化的采样率。具体采样过程代码实现和分析，请参见本文第8小节的介绍） 

到此我陷入了解这道题最大的困境，我的脱困思路是：既然保证K为整数的难点在于L必须是2的整数次幂，而L为2的整数次幂是由于快速算法FFT带来的。那不如干脆放弃快速傅里叶算法FFT，回到离散傅里叶变换DFT的本源，直接计算指定频点上的频域值（含幅度谱值和相位谱值）。另一方面，解本题也并不需要知道奈奎斯特频率以下所有谱线的信息，我只要计算C信号在10KHz、15KHz、20KHz……95KHz和100KHz等固定谱线上的幅度和相位即可！这也大大缓解了不使用FFT快速算法带来的计算量增加问题。下面是我所采用的具体算法，如有不妥之处欢迎读者和网友指正： 

(1)固定以系统可以达到的最高采样率1MSPS进行A/D和D/A转换，同时将算法使用的数据长度L固定为1000个采样点（即时长为1ms的信号）。

(2)在内存中生成采样率为1MSPS，长度L=1000个采样点，频率为10KHz、15KHz、20KHz……95KHz和100KHz的（共19种）正信号:

 

以及19种余弦信号： 

(3)当需要校正相位误差时（假设此时可称为t0），就同t0时刻开始对C信号实施长度为L=1000个采样点，历时1ms的采样。采样结果记为：  

 

(4)用采样结果分别与19种不同频率的正/余弦信号实施相关计算： 

(5)根据离散傅里叶变换（DFT）理论，由于{S_j,i}和{C_j,i}是正交的，上述两个计算结果 与 之比等于实际信号 中频率为 的信号的初始相位 的正切值。即可以通过下式计算初始相位：

(6)当然，根据题意，混合信号C中只有可能混有两种频率的信号。找到这两种频率的方法很简单，就是在计算（4）式之前，先通过功率谱计算公式： 

得到所有可能得19种信号的功率谱，再找到P_j中最大的两个所对应的频率即可。 

(7) 将P_j最大的两个频率的计算出来就是C中所含的A、B两个信号的相位差了。

当我把这个算法写完，发现这不就是软件无线电（Software Defined Radio(SDR）中的数字正交解调吗！仔细想想，好像还真是，这算法就是数字相位解调（Phase demodulation）——只是载波频率略低，以方便参赛者可以通过MCU验证软件无线电算法。这道题的命题老师应该是谙熟通信的专家，此题既帮助学生深刻的理解了DFT的定义，又有力的引出了软件无线电的算法核心思想，可谓用心良苦。

2）相位误差校正的实现

(1)相位同步A/D采集

实现相位误差校正首先要计算输入信号A与重建信号A’，以及B与B’之间的相位差（以下描述以A与A’为例介绍）。如前所述，计算A与A’相位差的方法是将长度为L个点（含有整数个信号周期）的A与A’信号代入公式（4）。由于A’是算法自己生成的信号，可以直接从生成信号的表格中读取，不需要进行A/D采集。为简便起见，从表格中读取参与（4）式计算的A’信号是从0相位开始的，因此对A信号的采样也需要在D/A重建的A’信号刚好播放到0相位时开始。所以图5所示的程序流程图右侧的循环相位校正部分的开始，就是循环等待到D/A输出的A’信号相位达到0相位的地方。举例来讲，我开通了DMA通道传输完成中断（DMA控制D/A转换的代码实现方式请参见本文第8小节的介绍），并在中断服务程序中对全局变量do_ad_flag1置位。中断服务程序代码如下所示：

void DMA2_Channel3_IRQHandler(void){    
    if(DMA_GetITStatus(DMA2_IT_TC3)) {//判断通道2是否传输完成
       DMA_ClearITPendingBit(DMA2_IT_TC3);    //清除通道2传输完成标志位
        do_ad_flag1 = 0x55;//把全局标志置位，表示可以同步开始AD转换了
     }
}

这样当主程序代码对do_ad_flag1清零后，只需要等待到下一次do_ad_flag1被置位的时刻，就是A’信号播放一轮L个点结束，也就是0相位的时刻。随后就可以启动DMA控制的A/D转换，实现长度为L的信号采样（DMA控制A/D转换的代码实现方式请参见本文第8小节的介绍）。

do_ad_flag1 = 0;//标志启动AD标志清零，等待完成一次DAC的DMA传输中断中置位
while(do_ad_flag1 == 0);//一直等待到一轮DA转换的DMA完成，才同步开启AD转换的DMA
ADC_DMA_Config(); //启动ADC的DMA

(2)离散余弦变换中相关（correlation）算法的实现

根据图5所示的程序整体流程图，公式（2）、（3）所示的离散余弦变换需要在程序执行过程中不断执行，以计算A和B信号的初始相位。因此执行公式（2）、（3）算法将占去CPU的最大计算带宽，需要尽量优化，以缩短执行两次相位校正的时间间隔，达到使输出的重建信号和输入信号尽量相位同步，避免示波器信号“显示漂移”的目的。好在A和B可选的频率只有19种，且可以在程序开始运行伊始就确定A和B的频率（对应图5中左侧的流程，其算法将在本文第6小节详细介绍）。对A’进行相位校正时只需要计算一次公式（2）、（3）即可，即（2）、（3）中的j值是固定的。我对公式（2）、（3）的实现代码如下所示：

double corr1000_200(short* data, short* mask) {
    double res = 0;
    short i;
    for(i = 0; i < 1000 ; i++)
        res = res + (data[i] - 2048) * (mask[i%200] - 2048);//注意模版要减去直流偏置
    return res;
}

其中，采集得到的信号T_i存储在指针*data所指向的数组中，固定频率的正、余弦信号SC_j和CC_j存储在*mask所指向的数组中。由图3所示的电路原理图可知，二者为了满足MCU的A/D和D/A的要求，都是加入了Vcc/2的直流偏置电压，在数字域对应的数值为2¹²÷2=2048，我在计算之前减去了这个直流偏置。另外，由于正、余弦信号表格中的数据SC_j和CC_j具有周期性，我只制作了长度为200个点的表格，其他的4/5长度的数据则通过周期沿拓的方式产生，具体代码如上所示。至于为什么制作长度为200个点的表格理由较复杂，推导和计算过程详见本文第7小节。 

调用上面的代码，实现公式（2）、（3）的代码如下所示： 

//计算离散傅里叶变换在单个频率点上的实部和虚部
temp_double_sin = corr1000_200(ADCConvertedValue, DAC_SIN + 200*input_frq_index[0]);
temp_double_cos = corr1000_200(ADCConvertedValue, DAC_COS + 200*input_frq_index[0]);

 其中，信号频率的编号存储在数组元素input_frq_index[0]中，编号的原则是10KHz信号对应的编号为0, 15KHz信号对应的编号为1, 20KHz信号对应的编号为2,……，100KHz信号对应的编号为18。则编号（frq_index）与真实频率（f）的换算关系为：

 而19种不同频率的正、余弦数据表格的数据都存放在指针DAC_SIN和DAC_COS所指向的大表格中，并依次排放（具体数据参见本文第7小节）。那么指向所需信号数据的首地址的指针就是DAC_SIN + 200*input_frq_index[0]了。

在上述算法中，都有意识的使用了double数据类型，目的一是为了防止计算中的溢出现象；而是为了简化后续三角函数计算中的的数据类型转换。 

(3)相位差计算 

使用公式（4）计算重建信号A’、B’和A、B之间的相位差，其代码如下： 

temp_double_angle = atan2(temp_double_cos, temp_double_sin);

 注：正切函数图像如下所示

 

图9 正切函数图像

 计算正切值就是已知上图中的x值（也就是角度）来计算y值，每个x都能对应唯一的y值。但本题需要反过来使用正切函数，已知y值（也就是正切值），计算x的角度值时，就会出现多个角度都满足正切值为y的情况。只使用中间一根曲线的来搜寻对应的x值显然是不合理的。 

回到实际的物理意义，两个信号相位差可能为[-π，π]中的任意角度，只使用中间一根曲线的x轴的范围则只可能在 [-π/2，π/2]之间，显然有问题。 

问题发生在上图的正切函数图像的y轴正切值——只有一个值，即正切计算三角关系的“对边”和“临边”的比，这个值没法反映“对边”和“临边”究竟哪个为正数，哪个为负数。只能得出角度在[-π/2，π/2]之间的结论，也就是相角向量只能在第一和第四相限；而要使计算得到的相角在[-π，π]之间，就必须知道“对边”和“临边”的符号。

具体来说，就是计算相位差时调用的反正切公式只能是C语言标准库中有两个参数的double atan2(double y, double x)，而不能是只有一个参数的double atan(double x)。 

(4)实时相位校正 

上面得到的相位差是弧度制的，而嵌入式程序只能控制D/A重建的信号进行相对延迟，因此程序首先要将相位差折算为合理的延迟时间。具体步骤如下：

如前所述，上面得到的相位差的值域在[-π/2，π/2]之间，既有可能是相位延迟，也有可能是相位提前。当重建信号相对于输入信号的相位提前时，相位差为负数，只需要将该负相位的绝对值折算为延迟时间即可。但当重建信号的相位落后于输入信号时，相位差为正数，无法通过简单的直接延迟来提前重建信号的相位。一种可行的办法是利用重建信号的周期性，让其延迟2π减去相位差，即进行一个周期的求补操作。延迟相位折算的代码如下。 

double pos_angle;//转换为延迟0-2pi，全部为整数的角度
double temp_double;
short temp_short;
if(angle < 0)//atan2函数计算得到的相位差在-pi到pi，但当其小于0时，DAC输出的信号无法提前相位，只能将所有相位全部延迟一个周期，即进行一个周期的求补操作
        pos_angle = -angle;
else
        pos_angle = 2*3.1415926 - angle;

 其中的变量pos_angle的单位是弧度，将其折算成延迟时间的公式如（7）所示。其中的f为信号频率，frq_index为信号频率编号，其换算关系使用了公式（6）。

 为了实现精确相位调整，我直接对控制DAC和DMA重建波形的TIM6/TIM7进行延时（关于基于DMA的DAC信号重建的代码，请参见本文第8小节）。为获得最高的时间分辨率，TIM6/TIM7的时钟被配置可用的最高频率——系统时钟56MHz。（7）式的计算结果delay_time单位为秒（S），而实际使用的值却是delay_time中含有多少个TIM6/TIM7时钟周期，因此（7）式的结果应再乘以56M。化简后得到代码如下：

temp_double = 5600*pos_angle / (3.1415926*((double)frq_index + 2));
temp_short = temp_double + 0.5;//转换为整数的延迟时间，避免舍入误差
return temp_short;

 注：此处提醒读者注意，为在STM32F1系列上获得1MSPS的A/D采样率，本文所述的系统使用了56MHz的系统时钟。如果直接使用本文代码,一定要首先将系统时钟修改为56MHz！如果你使用的是意法半导体官方的固件库和程序模板，则需要再system_stm32f10x.c中，将系统时钟配置代码修改如下： 

/*!< Uncomment the line corresponding to the desired System clock (SYSCLK) 

   frequency (after reset the HSI is used as SYSCLK source) */ 

//#define SYSCLK_FREQ_HSE    HSE_Value 

//#define SYSCLK_FREQ_20MHz  20000000 

//#define SYSCLK_FREQ_36MHz  36000000

//#define SYSCLK_FREQ_48MHz  48000000

#define SYSCLK_FREQ_56MHz  56000000

//#define SYSCLK_FREQ_72MHz  72000000 

言归正传，上面返回的计算结果temp_short是整数型，这样做的目的是方便将其写入TIM6的定时寄存器。而在进行数据类型转换之前在temp_double基础上加了0.5是为了防止舍入误差，以实现浮点型向整数型转换过程中的四舍五入操作。

TIM6是DAC1重建A’的定时器，这里使用了一个对定时器定时时间延迟更精确的小技巧：直接向当前定时数值中增加需要延迟的时钟数其代码如下。（注：上一段代码中计算得到的延迟时间，此时被转存到了数组delay_times[]的各个元素中。）

temp_int = TIM6->CNT;
temp_int = temp_int + 65535 - delay_times[0];
TIM6->CNT = temp_int;

相对应的，如果要实现对B’信号的相位校正，就要向控制DAC2的采样间隔的TIM7中增加延迟时钟数了，这里不再给出源码。

本环节揭示的信号处理和嵌入式系统教学中的问题： 

a) 《信号与系统》和《数字信号处理》教学中对于相位谱的讲解不如功率谱深入，这可能是由于相位谱的物理意义不清晰造成的，但至少在通信领域相位谱具有非常明确的物理意义和使用价值，值得《信号与系统》和《数字信号处理》的任课教师考虑花费学时强调一下。教学中，本例可以作为很好的一个相位谱应用的实例。

b) 可以发现，最近几届电赛题正在回归电子“设计”本身。作为电子电工教学的“指挥棒”，赛题变得更加贴近我国本科教学的实际，回归基本原理、基本技能；放弃对新型元件、新方法的一味偏执追求。以本赛题为例：用离散傅里叶变换实现相位调制的解调是属于《信号与系统》、《数字信号处理》和《通信原理》等课程中需要学生理解和掌握的基本理论，对这些理论的掌握和应用能力就是求解此题的必备技能。但为了防止参赛者使用现成软件无线电方案（多基于FPGA或专用ASIC芯片），造成不公平现象，本赛题使用了较低的信号频率。这种做法有效的避免了参赛队将精力放在“术”的层面；而是把精力集中在学习本科阶段必须掌握的基本技能，以及提高解决复杂工程问题的能力，这些“道”层面的东西上。
这也提示我们作为基层电子电工教学工作者，不应该把过多的精力放在讲解FPGA使用的操作流程或者某个具体芯片的配置流程这些花哨的“雕虫小技”，而应该不断抓牢基本原理、基本技能和解决复杂工程问题这些重点上。

c) 学生缺乏对算法复杂度的估计能力，这一方面源于学生的工程不足。但更重要的原因是“双师型”教师缺乏，没有相当的工程经验的教师很难对自己设计算法的计算量进行有效估计，更难交给学生对算法复杂度进行估算的习惯和方法。

6.  输入信号C中混合的A、B的频率和种类的判断

 1）频率的判断

C中混杂的两个分量A、B的频率只有可能是10KHz、15KHz……100KHz等共19种，所以对大多数参赛者而言，识别它们都不是一件特别困难的事。有同学使用FFT，当然也可以使用相关计算DFT。我自己在做的时候犹豫了一下，后来还是胡广书老师书中的话（图8）影响了我——如果使用长度L为2的整数次幂（如1024）点的FFT，则L内的数据对于所有19种频率都不刚好含有整数个周期，那么势必造成转换结果的频谱能量泄露。也就是不能以功率谱中单个点的幅度来代表信号在某个频率上的总能量，而要使用频谱中至少2-3个相邻点的功率之和，才能代表输入信号在单个频率上的功率。这种算法对于只有正弦输入的“基本要求”部分没有多大影响。但如果考虑“发挥部分”可能的三角波输入，就不得不考虑三角波在奇次谐波上的谐波失真和基频信号的叠加可能对功率谱的影响。

因此我决定仍然采用L=1000的DFT来计算信号在19个指定频率频点上的功率。此时由于时域长度为L=1000的信号中，含有的所有信号都是整数个周期，信号的所有频谱只会有一根谱线。在计算时无需像使用FFT算法一样，取相邻几个频点的功率，得到的功率值将更准确；在输入信号有可能是三角波时也不容易和其三次谐波混淆。加上前面所述的相位差算法中已经有了实现（2）、（3）两式的代码，只要再增加实现（5）式的代码即可。

//计算IQ分量的模的长度
double modulus(double I, double Q){
    double temp_double;
    temp_double = sqrt(I*I + Q*Q);
    return temp_double;
}

DFT计算两个输入信号频率的代码罗列与下面，供网友参考。

void cal_2frqs(unsigned char* frq)//计算参与叠加的两个正弦信号的频率
{
    double corr_sin_double;//用于存放C与所有19种正弦相关计算的结果
    double corr_cos_double;//用于存放C与所有19种余弦相关计算的结果
    double mag[19];//用于存放输入信号与所有19种信号相关计算的模
    double temp_double[19],sort_temp;
    unsigned char i,j,temp_char;
    for(i = 0; i < 19 ; i++)
    {
        corr_sin_double = corr1000_200(ADCConvertedValue, DAC_SIN + 200*i);
        corr_cos_double = corr1000_200(ADCConvertedValue, DAC_COS + 200*i);
        mag[i] = modulus(corr_sin_double, corr_cos_double);
        temp_double[i] = mag[i]; 
    }
    //冒泡法对大小进行排序
    for (i = 0; i < 19 - 1; i++)
    {
        for (j = 0; j < 19 - i - 1; j++){
                if (temp_double[j] < temp_double[j+1]){
                        sort_temp = temp_double[j];
                        temp_double[j] = temp_double[j+1];
                        temp_double[j+1] = sort_temp;
                }
            }
    }
    //寻找最大值序号
    for (i = 0; i < 19; i++){
        if(mag[i] == temp_double[0])//找到最大数
            frq[0] = i;
        if(mag[i] == temp_double[1])//找到第二大数
            frq[1] = i;
    }
    if(frq[0] > frq[1]) {//把频率低的放在第一个位置，如果颠倒就换一下
        temp_char = frq[0];
        frq[0] = frq[1];
        frq[1] = temp_char;
    }
}

其中，值得注意的是计算出所有的19个相关结果的模后，我用冒泡法对这些模进行了排序。而这个冒泡法排序的代码是百度“文心一言”写的，不得不说这提高了对数据结构不太熟悉的人的编程速度。然后逐一搜索找到了最大的两个的模所在的频率编号i，并把它们放到了数组frq[ ]中。啰嗦一句，至此放在frq[ ]中的编号和实际频率之间的对应关系为公式（6）。而上述函数真实通过向形式参数指针所指向的数组写入数值的方式向外传递计算结果的。

该函数的最后，还对frq[0]和frq[1]中存储的频率编号进行了调整，使较低频率的信号放在frq[0]中，这样做的目的使D/A转换器1输出的信号频率永远小于D/A转换器2输出的信号频率。

2）信号种类（正弦或三角波）的判断

“发挥部分”要求参赛者能够判断A和B信号究竟是正弦波还是三角波。从理论上讲这也并不困难——三角波在及次谐波上存在谐波分量。为简便起见，我选择对3次谐波的能量进行计算——并将其与基频功率进行比较，当达到一定比例时就认为这些频点上存在谐波失真，也就是输入是三角波。这就要求程序能够计算频率在30KHz，45KHz，……，250KHz，275KHz，300KHz频率位置的功率谱，如果继续使用前述的不带快速算法的DFT，就还要给出上面这些频点的正、余弦数据，也就需要在现有的正、余弦数据表的基础上，还要增加100KHz频率以上的信号数据表。

而使用FFT则可以直接得到奈奎斯特频率以下的所有谱线功率，因此我偷了个懒，直接使用了FFT。为使信号长度为2的整数次幂，以方便FFT计算，我将DMA控制A/D采样的长度增加到1024。由于1024点的数据不可能含有整数个波形，频谱的能量也不会集中在一根谱线上。下表所示的是FFT长度为1024个点时，我们感兴趣的三次谐波谱线在FFT结果数组中的位置。

 

表1 三次谐波频率在FFT结果中的位置

 

三次谐波频率（KHz）	三次谐波频率是奈奎斯特频率的倍数	1024点FFT所在的点（位置）
30	0.06	30.72
45	0.09	46.08
60	0.12	61.44
75	0.15	76.8
90	0.18	92.16
105	0.21	107.52
120	0.24	122.88
135	0.27	138.24
150	0.3	153.6
165	0.33	168.96
180	0.36	184.32
195	0.39	199.68
210	0.42	215.04
225	0.45	230.4
240	0.48	245.76
255	0.51	261.12
270	0.54	276.48
285	0.57	291.84
300	0.6	307.2
	频率分辨率（Hz）	0.9765625

可以发现，感兴趣的频点并不一定在整数点上。为防止相邻谱线中的能量泄露导致的问题，我对感兴趣频率谱线位置周边3个相邻点的功率进行了求和。

不过上述算法存在一个重大漏洞：当被测信号C中混合的A、B两个单频信号中的一个的频率刚好是另一个信号频率的3倍时（例如A为10KHz，B为30KHz时），低频信号的3次谐波会和较高频率信号的能量混淆在一起，从而使算法是无法判断低频信号是三角波还是正弦波的。但幸运的是，根据题意这种情况只会出现在三次谐波频率低于最高频率100KHz的情况下。具体来说会在（1）低频信号是10KHz，高频信号是30KHz；（2）低频信号是15KHz，高频信号是45KHz；（3）低频信号是20KHz，高频信号是60KHz；（4）低频信号是25KHz，高频信号是75KHz；（5）低频信号是30KHz，高频信号是90KHz。共五种情况下出现这个问题。

我在上述五种情况下，对低频信号的五次谐波进一步进