前言
由于本人接触的项目，DSP板没有留通信接口给PC端，因此要对系统做辨识，获取数据的唯一渠道就是直接把DSP内存里的数据通过CCS拷贝出来。本文重点讨论在DSP内存有限的情况下，高效地完成系统辨识的方法。
在做这部分工作前，需要学习两个前置知识，即
DSP_CCS7实现变量的导出与MatLAB读取_江湖上都叫我秋博的博客-CSDN博客
DSP_定义一个大的全局数组_探索之路_江湖上都叫我秋博的博客-CSDN博客
这两部分内容是必须要懂的，不然后续的内容，容易卡壳。
正文
低频大数据问题
假设系统采样率为5kHz，要辨识0.1Hz处系统特性，需要存储5wx2 = 10w个点（5w输入序列，5w输出序列）。前面文章DSP_定义一个大的全局数组_探索之路_江湖上都叫我秋博的博客-CSDN博客中已经探讨了，28377D最多存2.4w左右个点，显然是不够的。因此我们选择在低频做降采样率处理，用1kHz的采样率，0.1Hz需要1wx2=2w个点，刚好合适。 而更高的频率比如1Hz，我们则继续用5kHz的采样率。
频率点简化及高频数据合并
对于频段0.1~1000Hz常规的系统辨识一般来200个频率点都不过分， 那难到我们要去点CCS点200次，存200次数据？ 显然这样不够高效，因此对于高频部分，我们把多个频率点的数据都合并存到大小为1w的数组中（输入序列一个、输出序列一个）。另外，经过实验，只有能被采样率整除的频率点，它的数据点数才是整数，该频率点下的辨识更加准确。因此经过层层筛选， 200个频率被我筛到39个频率点。 基于上述的两个简化思路。 用MatLAB设计了分配策略。代码如下：
data_allocation.m
clear all;
close all;
clc;
% 频率 采样率 点数 周期数
FreqArray = [
0.1, 1000, 10000, 2;
0.125, 1000, 8000, 2;
0.16, 1000, 6250, 2;
0.2, 1000, 5000, 2;
0.25, 1000, 4000, 2;
0.3125, 1000, 3200, 2;
0.4, 1000, 2500, 2;
0.5, 1000, 2000, 2;
0.625, 1000, 1600, 4;
0.8, 1000, 1250, 4;
1, 5000, 5000, 4;
1.25, 5000, 4000, 4;
1.5625, 5000, 3200, 4;
2, 5000, 2500, 4;
2.5, 5000, 2000, 4;
3.125, 5000, 1600, 6;
4, 5000, 1250, 6;
5, 5000, 1000, 6;
6.25, 5000, 800, 6;
7.8125, 5000, 640, 6;
10, 5000, 500, 10;
12.5, 5000, 400, 10;
15.625, 5000, 320, 10;
20, 5000, 250, 10;
25, 5000, 200, 10;
31.25, 5000, 160, 10;
40, 5000, 125, 10;
50, 5000, 100, 10;
62.5, 5000, 80, 10;
78.125, 5000, 64, 10;
100, 5000, 50, 20;
125, 5000, 40, 20;
156.25, 5000, 32, 20;
200, 5000, 25, 20;
250, 5000, 20, 20;
312.5, 5000, 16, 20;
500, 5000, 10, 40;
625, 5000, 8, 40;
1000, 5000, 5, 40;
];
FreqNum = length(FreqArray(:,1));
% 数据校验
for i=1 : 1 : FreqNum
if (FreqArray(i, 2) / FreqArray(i, 1)) ~= FreqArray(i, 3)
fprintf('%f is wrong!\r\n',FreqArray(i, 1));
else
fprintf('%f is ok!\r\n',FreqArray(i, 1));
end
end
%-------------------------------------------------------------------------------------
% 频响数据的排列分配
% 数组总长度
DSP_Array_Length = 10000;
Current_Length = 0;
AllocationResult = zeros(FreqNum,4);
% 第一列：频率 第二列：频响数据序号
% 第三列：起始位置 第四列：截至位置
DataIndex = 1;
for i=1 : 1 : FreqNum
if (Current_Length + FreqArray(i, 3)) <= DSP_Array_Length
StartPosition = Current_Length + 1;
Current_Length = Current_Length + FreqArray(i, 3);
StopPosition = Current_Length;
AllocationResult(i,1) = FreqArray(i, 1);
AllocationResult(i,2) = DataIndex;
AllocationResult(i,3) = StartPosition;
AllocationResult(i,4) = StopPosition;
elseif FreqArray(i, 3) <= DSP_Array_Length
DataIndex = DataIndex + 1;
StartPosition = 1;
Current_Length = FreqArray(i, 3);
StopPosition = Current_Length;
AllocationResult(i,1) = FreqArray(i, 1);
AllocationResult(i,2) = DataIndex;
AllocationResult(i,3) = StartPosition;
AllocationResult(i,4) = StopPosition;
else
fprintf('频率点：%f的数据点数超过了DSP的数组总长度 \r\n',FreqArray(i, 1));
end
end
%-------------------------------------------------------------------------------------
%--------------------------------------------------------------------------------------
% DSP代码生成
para_path = 'AutoGeneratedCode.txt';
fp = fopen(para_path, 'w');
fprintf(fp,'**请把下面的代码拷贝到dsp工程项目的sys_identify_init函数中**\n');
Current_Length = 0;
DataIndex = 1;
for i=1 : 1 : FreqNum
if (Current_Length + FreqArray(i, 3)) <= DSP_Array_Length
StartPosition = Current_Length + 1;
Current_Length = Current_Length + FreqArray(i, 3);
StopPosition = Current_Length;
AllocationResult(i,1) = FreqArray(i, 1);
AllocationResult(i,2) = DataIndex;
AllocationResult(i,3) = StartPosition;
AllocationResult(i,4) = StopPosition;
fprintf(fp,' FreqPoints[%d].freq = %f;\n', i-1, FreqArray(i, 1));
fprintf(fp,' FreqPoints[%d].sampling_rate = %d;\n', i-1, FreqArray(i, 2));
fprintf(fp,' FreqPoints[%d].data_number = %d;\n', i-1, FreqArray(i, 3));
fprintf(fp,' FreqPoints[%d].data_index = %d;\n', i-1, DataIndex);
fprintf(fp,' FreqPoints[%d].start_position = %d;\n', i-1, StartPosition - 1);
fprintf(fp,' FreqPoints[%d].stop_position = %d;\n', i-1, StopPosition - 1);
fprintf(fp,' FreqPoints[%d].period_number = %d;\n', i-1, FreqArray(i, 4));
fprintf(fp,' FreqPoints[%d].freq_div_addr = &iFreq%dHz;\n\n',i-1,FreqArray(i, 2));
elseif FreqArray(i, 3) <= DSP_Array_Length
DataIndex = DataIndex + 1;
StartPosition = 1;
Current_Length = FreqArray(i, 3);
StopPosition = Current_Length;
AllocationResult(i,1) = FreqArray(i, 1);
AllocationResult(i,2) = DataIndex;
AllocationResult(i,3) = StartPosition;
AllocationResult(i,4) = StopPosition;
fprintf(fp,' FreqPoints[%d].freq = %f;\n', i-1, FreqArray(i, 1));
fprintf(fp,' FreqPoints[%d].sampling_rate = %d;\n', i-1, FreqArray(i, 2));
fprintf(fp,' FreqPoints[%d].data_number = %d;\n', i-1, FreqArray(i, 3));
fprintf(fp,' FreqPoints[%d].data_index = %d;\n', i-1, DataIndex);
fprintf(fp,' FreqPoints[%d].start_position = %d;\n', i-1, StartPosition - 1);
fprintf(fp,' FreqPoints[%d].stop_position = %d;\n', i-1, StopPosition - 1);
fprintf(fp,' FreqPoints[%d].period_number = %d;\n', i-1, FreqArray(i, 4));
fprintf(fp,' FreqPoints[%d].freq_div_addr = &iFreq%dHz;\n\n',i-1,FreqArray(i, 2));
else
fprintf('频率点：%f的数据点数超过了DSP的数组总长度 \r\n',FreqArray(i, 1));
end
end
fclose(fp);
open(para_path);
根据代码运行结果，系统实现一次完整的辨识，只需要拷贝8次结果。 相比于200次拷贝优化程度已经较高了。 另外，代码中还包含了dsp部分代码的生成，这部分代码可以拷贝到dsp的项目当中，这对开发人员（我）是福音啊。
DSP相关代码
要完成系统辨识，dsp这边的代码肯定是必不可少的。
1、针对不同的采样率，我们需要用定时器（CPU Timer）进行分频实现。
2、针对每个频率点需要用到信息，拷贝序列的安排，拷贝上述.m文件生成的代码即可完成初始化
3、具体的整体实现流程也并不复杂，详情看代码就行了，自行体会。
main.c
/*
* main.c
*
* Created on: 2022年12月8日
* Author: YFR-HQN
*/
#include <main.h>
int16 iFreqCnt = 0;
int16 iFreq1Hz = 0;
int16 iFreq2Hz = 0;
int16 iFreq20Hz = 0;
int16 iFreq40Hz = 0;
int16 iFreq50Hz = 0;
int16 iFreq100Hz = 0;
int16 iFreq200Hz = 0;
int16 iFreq250Hz = 0;
int16 iFreq500Hz = 0;
int16 iFreq1000Hz = 0;
int16 iFreq2000Hz = 0;
int16 iFreq2500Hz = 0;
int16 iFreq5000Hz = 0;
void main(void)
{
default_init();
interrupt_init();
CpuTimer0Regs.TCR.bit.TSS = 0; // 使能CPUTimer0
runtime_init();
sys_identify_init();
sys_idenfify();
}
void default_init(void){
InitSysCtrl();
DINT;
InitPieCtrl();
InitGpio();
IER = 0x0000;
IFR = 0x0000;
InitPieVectTable();
#ifdef _STANDALONE
#ifdef _FLASH
IPCBootCPU2(C1C2_BROM_BOOTMODE_BOOT_FROM_FLASH);
#else
IPCBootCPU2(C1C2_BROM_BOOTMODE_BOOT_FROM_RAM);
#endif
#endif
InitCpuTimers(); // Initialize CPU Timers
}
void interrupt_init(void){
// CPU Timer0
ConfigCpuTimer(&CpuTimer0, 200, 1000000/SYS_BASIC_FREQUENCY);
CpuTimer0Regs.TCR.bit.TIE = 1; // 1: The CPU-Timer interrupt is enabled
// TIMER0 INT1.7 [Page102: Table 3-2.PIE Channel Mapping]
EALLOW;
PieVectTable.TIMER0_INT = &CpuTimer0ISR; // Specify the interrupt service routine
EDIS;
IER |= M_INT1; // Enable CPU Level interrupt
PieCtrlRegs.PIEIER1.bit.INTx7 = 1; // Enable PIE Level interrupt
PieCtrlRegs.PIECTRL.bit.ENPIE = 1;
EINT; // 使能全局中断
ERTM; // 使能实时仿真中断
}
// Generate frequency division signal
interrupt void CpuTimer0ISR(void){
DivideFreq();
PieCtrlRegs.PIEACK.all = PIEACK_GROUP1;
}
void DivideFreq(void){
if(iFreqCnt >= CNT1HZ){
iFreq1Hz = 1;
iFreqCnt = 0;
}else{
iFreq1Hz = 0;
}
if((iFreqCnt%CNT2HZ)==0){
iFreq2Hz = 1;
}else{
iFreq2Hz = 0;
}
if((iFreqCnt%CNT20HZ)==0){
iFreq20Hz = 1;
}else{
iFreq20Hz = 0;
}
if((iFreqCnt%CNT40HZ)==0){
iFreq40Hz = 1;
}else{
iFreq40Hz = 0;
}
if((iFreqCnt%CNT50HZ)==0){
iFreq50Hz = 1;
}else{
iFreq50Hz = 0;
}
if((iFreqCnt%CNT100HZ)==0){
iFreq100Hz = 1;
}else{
iFreq100Hz = 0;
}
if((iFreqCnt%CNT200HZ)==0){
iFreq200Hz = 1;
}else{
iFreq200Hz = 0;
}
if((iFreqCnt%CNT250HZ)==0){
iFreq250Hz = 1;
}else{
iFreq250Hz = 0;
}
if((iFreqCnt%CNT500HZ)==0){
iFreq500Hz = 1;
}else{
iFreq500Hz = 0;
}
if((iFreqCnt%CNT1000HZ)==0){
iFreq1000Hz = 1;
}else{
iFreq1000Hz = 0;
}
/*
if((iFreqCnt%CNT2000HZ)==0){
iFreq2000Hz = 1;
}else{
iFreq2000Hz = 0;
}
*/
if((iFreqCnt%CNT2500HZ)==0){
iFreq2500Hz = 1;
}else{
iFreq2500Hz = 0;
}
if((iFreqCnt%CNT5000HZ)==0){
iFreq5000Hz = 1;
}else{
iFreq5000Hz = 0;
}
iFreqCnt++;
}
main.h
/*
* main.h
*
* Created on: 2022年12月8日
* Author: YFR-HQN
*/
#ifndef PROGRAM_MAIN_H_
#define PROGRAM_MAIN_H_
#include <F28x_Project.h>
#include <F2837xD_Ipc_drivers.h>
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <runtime.h>
#include <sys_identify.h>
#include <control.h>
void default_init(void);
void interrupt_init(void);
interrupt void CpuTimer0ISR(void);
void DivideFreq(void);
#define SYS_BASIC_FREQUENCY 5000
extern int16 iFreqCnt;
extern int16 iFreq1Hz;
extern int16 iFreq2Hz;
extern int16 iFreq20Hz;
extern int16 iFreq40Hz;
extern int16 iFreq50Hz;
extern int16 iFreq100Hz;
extern int16 iFreq200Hz;
extern int16 iFreq250Hz;
extern int16 iFreq500Hz;
extern int16 iFreq1000Hz;
extern int16 iFreq2500Hz;
extern int16 iFreq5000Hz;
#define CNT1HZ 5000
#define CNT2HZ 2500
#define CNT20HZ 250
#define CNT40HZ 125
#define CNT50HZ 100
#define CNT100HZ 50
#define CNT200HZ 25
#define CNT250HZ 20
#define CNT500HZ 10
#define CNT1000HZ 5
#define CNT2500HZ 2
#define CNT5000HZ 1
#endif /* PROGRAM_MAIN_H_ */
sys_identify.c
/*
* sys_identify.c
*
* Created on: 2023年2月6日
* Author: YFR-HQN
*/
#include <main.h>
float sys_identify_in[10000];
float sys_identify_out[10000];
// float sys_identify_in = 0.0;
int FreqPonitNum = 39;
struct FreqPoint FreqPoints[39];
void sys_idenfify(void){
int i = 0;
int j = 0;
int k = 0;
float t = 0.0;
float delta_t = 0;
float pi = 3.141592653589793;
// 频率点
for(i = 0; i < FreqPonitNum; i++){
delta_t = __divf32(1,(float)FreqPoints[i].sampling_rate);
// 测试周期
for(j = 0; j < FreqPoints[i].period_number; j++){
t = 0.0;
for(k = 0;k < FreqPoints[i].data_number;k++){
while(*FreqPoints[i].freq_div_addr != 1);
*FreqPoints[i].freq_div_addr = 0;
// sys_identify_in = __sin(2 * pi * FreqPoints[i].freq * t);
if(j == (FreqPoints[i].period_number-1)){
sys_identify_in[FreqPoints[i].start_position + k] = __sin(2 * pi * FreqPoints[i].freq * t);
sys_identify_out[FreqPoints[i].start_position + k] = 10 * __sin(2 * pi * FreqPoints[i].freq * t + pi/3.0);
}
t = t + delta_t;
}
}
if((FreqPoints[i].start_position + k - 1)!=FreqPoints[i].stop_position){
asm (" ESTOP0"); // 数据分配错误
}
if(i == (FreqPonitNum-1)){
asm (" ESTOP0"); // 断点，保存数据
}else if(FreqPoints[i+1].data_index != FreqPoints[i].data_index){
asm (" ESTOP0"); // 断点，保存数据
}
}
}
void sys_identify_init(void){
FreqPoints[0].freq = 0.100000;
FreqPoints[0].sampling_rate = 1000;
FreqPoints[0].data_number = 10000;
FreqPoints[0].data_index = 1;
FreqPoints[0].start_position = 0;
FreqPoints[0].stop_position = 9999;
FreqPoints[0].period_number = 2;
FreqPoints[0].freq_div_addr = &iFreq1000Hz;
FreqPoints[1].freq = 0.125000;
FreqPoints[1].sampling_rate = 1000;
FreqPoints[1].data_number = 8000;
FreqPoints[1].data_index = 2;
FreqPoints[1].start_position = 0;
FreqPoints[1].stop_position = 7999;
FreqPoints[1].period_number = 2;
FreqPoints[1].freq_div_addr = &iFreq1000Hz;
FreqPoints[2].freq = 0.160000;
FreqPoints[2].sampling_rate = 1000;
FreqPoints[2].data_number = 6250;
FreqPoints[2].data_index = 3;
FreqPoints[2].start_position = 0;
FreqPoints[2].stop_position = 6249;
FreqPoints[2].period_number = 2;
FreqPoints[2].freq_div_addr = &iFreq1000Hz;
FreqPoints[3].freq = 0.200000;
FreqPoints[3].sampling_rate = 1000;
FreqPoints[3].data_number = 5000;
FreqPoints[3].data_index = 4;
FreqPoints[3].start_position = 0;
FreqPoints[3].stop_position = 4999;
FreqPoints[3].period_number = 2;
FreqPoints[3].freq_div_addr = &iFreq1000Hz;
FreqPoints[4].freq = 0.250000;
FreqPoints[4].sampling_rate = 1000;
FreqPoints[4].data_number = 4000;
FreqPoints[4].data_index = 4;
FreqPoints[4].start_position = 5000;
FreqPoints[4].stop_position = 8999;
FreqPoints[4].period_number = 2;
FreqPoints[4].freq_div_addr = &iFreq1000Hz;
FreqPoints[5].freq = 0.312500;
FreqPoints[5].sampling_rate = 1000;
FreqPoints[5].data_number = 3200;
FreqPoints[5].data_index = 5;
FreqPoints[5].start_position = 0;
FreqPoints[5].stop_position = 3199;
FreqPoints[5].period_number = 2;
FreqPoints[5].freq_div_addr = &iFreq1000Hz;
FreqPoints[6].freq = 0.400000;
FreqPoints[6].sampling_rate = 1000;
FreqPoints[6].data_number = 2500;
FreqPoints[6].data_index = 5;
FreqPoints[6].start_position = 3200;
FreqPoints[6].stop_position = 5699;
FreqPoints[6].period_number = 2;
FreqPoints[6].freq_div_addr = &iFreq1000Hz;
FreqPoints[7].freq = 0.500000;
FreqPoints[7].sampling_rate = 1000;
FreqPoints[7].data_number = 2000;
FreqPoints[7].data_index = 5;
FreqPoints[7].start_position = 5700;
FreqPoints[7].stop_position = 7699;
FreqPoints[7].period_number = 2;
FreqPoints[7].freq_div_addr = &iFreq1000Hz;
FreqPoints[8].freq = 0.625000;
FreqPoints[8].sampling_rate = 1000;
FreqPoints[8].data_number = 1600;
FreqPoints[8].data_index = 5;
FreqPoints[8].start_position = 7700;
FreqPoints[8].stop_position = 9299;
FreqPoints[8].period_number = 4;
FreqPoints[8].freq_div_addr = &iFreq1000Hz;
FreqPoints[9].freq = 0.800000;
FreqPoints[9].sampling_rate = 1000;
FreqPoints[9].data_number = 1250;
FreqPoints[9].data_index = 6;
FreqPoints[9].start_position = 0;
FreqPoints[9].stop_position = 1249;
FreqPoints[9].period_number = 4;
FreqPoints[9].freq_div_addr = &iFreq1000Hz;
FreqPoints[10].freq = 1.000000;
FreqPoints[10].sampling_rate = 5000;
FreqPoints[10].data_number = 5000;
FreqPoints[10].data_index = 6;
FreqPoints[10].start_position = 1250;
FreqPoints[10].stop_position = 6249;
FreqPoints[10].period_number = 4;
FreqPoints[10].freq_div_addr = &iFreq5000Hz;
FreqPoints[11].freq = 1.250000;
FreqPoints[11].sampling_rate = 5000;
FreqPoints[11].data_number = 4000;
FreqPoints[11].data_index = 7;
FreqPoints[11].start_position = 0;
FreqPoints[11].stop_position = 3999;
FreqPoints[11].period_number = 4;
FreqPoints[11].freq_div_addr = &iFreq5000Hz;
FreqPoints[12].freq = 1.562500;
FreqPoints[12].sampling_rate = 5000;
FreqPoints[12].data_number = 3200;
FreqPoints[12].data_index = 7;
FreqPoints[12].start_position = 4000;
FreqPoints[12].stop_position = 7199;
FreqPoints[12].period_number = 4;
FreqPoints[12].freq_div_addr = &iFreq5000Hz;
FreqPoints[13].freq = 2.000000;
FreqPoints[13].sampling_rate = 5000;
FreqPoints[13].data_number = 2500;
FreqPoints[13].data_index = 7;
FreqPoints[13].start_position = 7200;
FreqPoints[13].stop_position = 9699;
FreqPoints[13].period_number = 4;
FreqPoints[13].freq_div_addr = &iFreq5000Hz;
FreqPoints[14].freq = 2.500000;
FreqPoints[14].sampling_rate = 5000;
FreqPoints[14].data_number = 2000;
FreqPoints[14].data_index = 8;
FreqPoints[14].start_position = 0;
FreqPoints[14].stop_position = 1999;
FreqPoints[14].period_number = 4;
FreqPoints[14].freq_div_addr = &iFreq5000Hz;
FreqPoints[15].freq = 3.125000;
FreqPoints[15].sampling_rate = 5000;
FreqPoints[15].data_number = 1600;
FreqPoints[15].data_index = 8;
FreqPoints[15].start_position = 2000;
FreqPoints[15].stop_position = 3599;
FreqPoints[15].period_number = 6;
FreqPoints[15].freq_div_addr = &iFreq5000Hz;
FreqPoints[16].freq = 4.000000;
FreqPoints[16].sampling_rate = 5000;
FreqPoints[16].data_number = 1250;
FreqPoints[16].data_index = 8;
FreqPoints[16].start_position = 3600;
FreqPoints[16].stop_position = 4849;
FreqPoints[16].period_number = 6;
FreqPoints[16].freq_div_addr = &iFreq5000Hz;
FreqPoints[17].freq = 5.000000;
FreqPoints[17].sampling_rate = 5000;
FreqPoints[17].data_number = 1000;
FreqPoints[17].data_index = 8;
FreqPoints[17].start_position = 4850;
FreqPoints[17].stop_position = 5849;
FreqPoints[17].period_number = 6;
FreqPoints[17].freq_div_addr = &iFreq5000Hz;
FreqPoints[18].freq = 6.250000;
FreqPoints[18].sampling_rate = 5000;
FreqPoints[18].data_number = 800;
FreqPoints[18].data_index = 8;
FreqPoints[18].start_position = 5850;
FreqPoints[18].stop_position = 6649;
FreqPoints[18].period_number = 6;
FreqPoints[18].freq_div_addr = &iFreq5000Hz;
FreqPoints[19].freq = 7.812500;
FreqPoints[19].sampling_rate = 5000;
FreqPoints[19].data_number = 640;
FreqPoints[19].data_index = 8;
FreqPoints[19].start_position = 6650;
FreqPoints[19].stop_position = 7289;
FreqPoints[19].period_number = 6;
FreqPoints[19].freq_div_addr = &iFreq5000Hz;
FreqPoints[20].freq = 10.000000;
FreqPoints[20].sampling_rate = 5000;
FreqPoints[20].data_number = 500;
FreqPoints[20].data_index = 8;
FreqPoints[20].start_position = 7290;
FreqPoints[20].stop_position = 7789;
FreqPoints[20].period_number = 10;
FreqPoints[20].freq_div_addr = &iFreq5000Hz;
FreqPoints[21].freq = 12.500000;
FreqPoints[21].sampling_rate = 5000;
FreqPoints[21].data_number = 400;
FreqPoints[21].data_index = 8;
FreqPoints[21].start_position = 7790;
FreqPoints[21].stop_position = 8189;
FreqPoints[21].period_number = 10;
FreqPoints[21].freq_div_addr = &iFreq5000Hz;
FreqPoints[22].freq = 15.625000;
FreqPoints[22].sampling_rate = 5000;
FreqPoints[22].data_number = 320;
FreqPoints[22].data_index = 8;
FreqPoints[22].start_position = 8190;
FreqPoints[22].stop_position = 8509;
FreqPoints[22].period_number = 10;
FreqPoints[22].freq_div_addr = &iFreq5000Hz;
FreqPoints[23].freq = 20.000000;
FreqPoints[23].sampling_rate = 5000;
FreqPoints[23].data_number = 250;
FreqPoints[23].data_index = 8;
FreqPoints[23].start_position = 8510;
FreqPoints[23].stop_position = 8759;
FreqPoints[23].period_number = 10;
FreqPoints[23].freq_div_addr = &iFreq5000Hz;
FreqPoints[24].freq = 25.000000;
FreqPoints[24].sampling_rate = 5000;
FreqPoints[24].data_number = 200;
FreqPoints[24].data_index = 8;
FreqPoints[24].start_position = 8760;
FreqPoints[24].stop_position = 8959;
FreqPoints[24].period_number = 10;
FreqPoints[24].freq_div_addr = &iFreq5000Hz;
FreqPoints[25].freq = 31.250000;
FreqPoints[25].sampling_rate = 5000;
FreqPoints[25].data_number = 160;
FreqPoints[25].data_index = 8;
FreqPoints[25].start_position = 8960;
FreqPoints[25].stop_position = 9119;
FreqPoints[25].period_number = 10;
FreqPoints[25].freq_div_addr = &iFreq5000Hz;
FreqPoints[26].freq = 40.000000;
FreqPoints[26].sampling_rate = 5000;
FreqPoints[26].data_number = 125;
FreqPoints[26].data_index = 8;
FreqPoints[26].start_position = 9120;
FreqPoints[26].stop_position = 9244;
FreqPoints[26].period_number = 10;
FreqPoints[26].freq_div_addr = &iFreq5000Hz;
FreqPoints[27].freq = 50.000000;
FreqPoints[27].sampling_rate = 5000;
FreqPoints[27].data_number = 100;
FreqPoints[27].data_index = 8;
FreqPoints[27].start_position = 9245;
FreqPoints[27].stop_position = 9344;
FreqPoints[27].period_number = 10;
FreqPoints[27].freq_div_addr = &iFreq5000Hz;
FreqPoints[28].freq = 62.500000;
FreqPoints[28].sampling_rate = 5000;
FreqPoints[28].data_number = 80;
FreqPoints[28].data_index = 8;
FreqPoints[28].start_position = 9345;
FreqPoints[28].stop_position = 9424;
FreqPoints[28].period_number = 10;
FreqPoints[28].freq_div_addr = &iFreq5000Hz;
FreqPoints[29].freq = 78.125000;
FreqPoints[29].sampling_rate = 5000;
FreqPoints[29].data_number = 64;
FreqPoints[29].data_index = 8;
FreqPoints[29].start_position = 9425;
FreqPoints[29].stop_position = 9488;
FreqPoints[29].period_number = 10;
FreqPoints[29].freq_div_addr = &iFreq5000Hz;
FreqPoints[30].freq = 100.000000;
FreqPoints[30].sampling_rate = 5000;
FreqPoints[30].data_number = 50;
FreqPoints[30].data_index = 8;
FreqPoints[30].start_position = 9489;
FreqPoints[30].stop_position = 9538;
FreqPoints[30].period_number = 20;
FreqPoints[30].freq_div_addr = &iFreq5000Hz;
FreqPoints[31].freq = 125.000000;
FreqPoints[31].sampling_rate = 5000;
FreqPoints[31].data_number = 40;
FreqPoints[31].data_index = 8;
FreqPoints[31].start_position = 9539;
FreqPoints[31].stop_position = 9578;
FreqPoints[31].period_number = 20;
FreqPoints[31].freq_div_addr = &iFreq5000Hz;
FreqPoints[32].freq = 156.250000;
FreqPoints[32].sampling_rate = 5000;
FreqPoints[32].data_number = 32;
FreqPoints[32].data_index = 8;
FreqPoints[32].start_position = 9579;
FreqPoints[32].stop_position = 9610;
FreqPoints[32].period_number = 20;
FreqPoints[32].freq_div_addr = &iFreq5000Hz;
FreqPoints[33].freq = 200.000000;
FreqPoints[33].sampling_rate = 5000;
FreqPoints[33].data_number = 25;
FreqPoints[33].data_index = 8;
FreqPoints[33].start_position = 9611;
FreqPoints[33].stop_position = 9635;
FreqPoints[33].period_number = 20;
FreqPoints[33].freq_div_addr = &iFreq5000Hz;
FreqPoints[34].freq = 250.000000;
FreqPoints[34].sampling_rate = 5000;
FreqPoints[34].data_number = 20;
FreqPoints[34].data_index = 8;
FreqPoints[34].start_position = 9636;
FreqPoints[34].stop_position = 9655;
FreqPoints[34].period_number = 20;
FreqPoints[34].freq_div_addr = &iFreq5000Hz;
FreqPoints[35].freq = 312.500000;
FreqPoints[35].sampling_rate = 5000;
FreqPoints[35].data_number = 16;
FreqPoints[35].data_index = 8;
FreqPoints[35].start_position = 9656;
FreqPoints[35].stop_position = 9671;
FreqPoints[35].period_number = 20;
FreqPoints[35].freq_div_addr = &iFreq5000Hz;
FreqPoints[36].freq = 500.000000;
FreqPoints[36].sampling_rate = 5000;
FreqPoints[36].data_number = 10;
FreqPoints[36].data_index = 8;
FreqPoints[36].start_position = 9672;
FreqPoints[36].stop_position = 9681;
FreqPoints[36].period_number = 40;
FreqPoints[36].freq_div_addr = &iFreq5000Hz;
FreqPoints[37].freq = 625.000000;
FreqPoints[37].sampling_rate = 5000;
FreqPoints[37].data_number = 8;
FreqPoints[37].data_index = 8;
FreqPoints[37].start_position = 9682;
FreqPoints[37].stop_position = 9689;
FreqPoints[37].period_number = 40;
FreqPoints[37].freq_div_addr = &iFreq5000Hz;
FreqPoints[38].freq = 1000.000000;
FreqPoints[38].sampling_rate = 5000;
FreqPoints[38].data_number = 5;
FreqPoints[38].data_index = 8;
FreqPoints[38].start_position = 9690;
FreqPoints[38].stop_position = 9694;
FreqPoints[38].period_number = 40;
FreqPoints[38].freq_div_addr = &iFreq5000Hz;
}
sys_identify.h
/*
* sys_identify.h
*
* Created on: 2023年2月6日
* Author: YFR-HQN
*/
#ifndef PROGRAM_SYSIDENTIFY_SYS_IDENTIFY_H_
#define PROGRAM_SYSIDENTIFY_SYS_IDENTIFY_H_
extern float sys_identify_out[10000];
extern float sys_identify_in[10000];
extern int FreqPonitNum;
extern struct FreqPoint FreqPoints[39];
struct FreqPoint{
float freq;
int sampling_rate;
int data_number;
int data_index;
int start_position;
int stop_position;
int period_number;
int16 * freq_div_addr;
};
void sys_identify_init(void);
void sys_idenfify(void);
#endif /* PROGRAM_SYSIDENTIFY_SYS_IDENTIFY_H_ */
注意：代码运行起来后，会在下面这部分代码进入自动进入中断
if(i == (FreqPonitNum-1)){
asm (" ESTOP0"); // 断点，保存数据
}else if(FreqPoints[i+1].data_index != FreqPoints[i].data_index){
asm (" ESTOP0"); // 断点，保存数据
}这时候，我们需要把分别把sys_identify_in数组和sys_identify_out数组导出来。
Tips: 在存数据的时候，把FreqPoints[i].data_index放到变量观察栏，根据这个变量判断当前的序列号是多少。 比如FreqPoints[i].data_index=1，那么就把sys_identify_in数组存为“../DspData/in1.dat"，把sys_identify_out存为“../DspData/out1.dat"，以此类推，..., FreqPoints[i].data_index=8时，sys_identify_in数组存为“../DspData/in8.dat"，sys_identify_out存为“../DspData/Out8.dat"。 在不停的存储数据时，除了修改文件名以外，还需要注意start_address那里需要不停的切换sys_identify_in/sys_identify_out。 看看下面截图。
8组输入输出序列存完后，如下所示：
MatLAB数据解析，系统辨识
得到了数据之后，就需要把上述的8组输入输出序列进行解析，得到其频域上的特性。具体的代码如下：
clear all; close all; clc;
FreqArray = [
0.1, 1000, 10000, 2;
0.125, 1000, 8000, 2;
0.16, 1000, 6250, 2;
0.2, 1000, 5000, 2;
0.25, 1000, 4000, 2;
0.3125, 1000, 3200, 2;
0.4, 1000, 2500, 2;
0.5, 1000, 2000, 2;
0.625, 1000, 1600, 4;
0.8, 1000, 1250, 4;
1, 5000, 5000, 4;
1.25, 5000, 4000, 4;
1.5625, 5000, 3200, 4;
2, 5000, 2500, 4;
2.5, 5000, 2000, 4;
3.125, 5000, 1600, 6;
4, 5000, 1250, 6;
5, 5000, 1000, 6;
6.25, 5000, 800, 6;
7.8125, 5000, 640, 6;
10, 5000, 500, 10;
12.5, 5000, 400, 10;
15.625, 5000, 320, 10;
20, 5000, 250, 10;
25, 5000, 200, 10;
31.25, 5000, 160, 10;
40, 5000, 125, 10;
50, 5000, 100, 10;
62.5, 5000, 80, 10;
78.125, 5000, 64, 10;
100, 5000, 50, 20;
125, 5000, 40, 20;
156.25, 5000, 32, 20;
200, 5000, 25, 20;
250, 5000, 20, 20;
312.5, 5000, 16, 20;
500, 5000, 10, 40;
625, 5000, 8, 40;
1000, 5000, 5, 40;
];
FreqNum = length(FreqArray(:,1));
for i=1 : 1 : FreqNum
if (FreqArray(i, 2) / FreqArray(i, 1)) ~= FreqArray(i, 3)
%fprintf('%f is wrong!\r\n',FreqArray(i, 1));
else
%fprintf('%f is ok!\r\n',FreqArray(i, 1));
end
end
% 数组总长度
DSP_Array_Length = 10000;
Current_Length = 0;
AllocationResult = zeros(FreqNum,4);
% 第一列表示频率 % 第二列表示频响数据序号
% 第三列起始位置 % 第四列截至位置
DataIndex = 1;
for i=1 : 1 : FreqNum
if (Current_Length + FreqArray(i, 3)) <= DSP_Array_Length
StartPosition = Current_Length + 1;
Current_Length = Current_Length + FreqArray(i, 3);
StopPosition = Current_Length;
AllocationResult(i,1) = FreqArray(i, 1);
AllocationResult(i,2) = DataIndex;
AllocationResult(i,3) = StartPosition;
AllocationResult(i,4) = StopPosition;
elseif FreqArray(i, 3) <= DSP_Array_Length
DataIndex = DataIndex + 1;
StartPosition = 1;
Current_Length = FreqArray(i, 3);
StopPosition = Current_Length;
AllocationResult(i,1) = FreqArray(i, 1);
AllocationResult(i,2) = DataIndex;
AllocationResult(i,3) = StartPosition;
AllocationResult(i,4) = StopPosition;
else
fprintf('频率点：%f的数据点数超过了DSP的数组总长度 \r\n',FreqArray(i, 1));
end
end
%-----------------------------------------------------------------------------------
%-----------------------------------------------------------------------------------
% 注意文件头的字节个数
%-----------------------------------------------------------------------------------
%-----------------------------------------------------------------------------------
fid = fopen('DspData/in1.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y1 = x{1};
fclose(fid);
fid = fopen('DspData/in2.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y2 = x{1};
fclose(fid);
fid = fopen('DspData/in3.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y3 = x{1};
fclose(fid);
fid = fopen('DspData/in4.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y4 = x{1};
fclose(fid);
fid = fopen('DspData/in5.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y5 = x{1};
fclose(fid);
fid = fopen('DspData/in6.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y6 = x{1};
fclose(fid);
fid = fopen('DspData/in7.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y7 = x{1};
fclose(fid);
fid = fopen('DspData/in8.dat');
fseek(fid, 21,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y8 = x{1};
fclose(fid);
y1 = y1';y2 = y2';y3=y3';y4=y4';y5=y5';y6=y6';y7=y7';y8=y8';
sys_identify_in = [y1;y2;y3;y4;y5;y6;y7;y8];
%-----------------------------------------------------------------------------------
%-----------------------------------------------------------------------------------
% 注意文件头的字节个数
%-----------------------------------------------------------------------------------
%-----------------------------------------------------------------------------------
fid = fopen('DspData/out1.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y1 = x{1};
fclose(fid);
fid = fopen('DspData/out2.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y2 = x{1};
fclose(fid);
fid = fopen('DspData/out3.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y3 = x{1};
fclose(fid);
fid = fopen('DspData/out4.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y4 = x{1};
fclose(fid);
fid = fopen('DspData/out5.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y5 = x{1};
fclose(fid);
fid = fopen('DspData/out6.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y6 = x{1};
fclose(fid);
fid = fopen('DspData/out7.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y7 = x{1};
fclose(fid);
fid = fopen('DspData/out8.dat');
fseek(fid, 23,-1); % 去文件头, 21字节，如 1651 1 80000006 0 100
x = textscan(fid, '%f');
y8 = x{1};
fclose(fid);
y1 = y1';y2 = y2';y3=y3';y4=y4';y5=y5';y6=y6';y7=y7';y8=y8';
sys_identify_out = [y1;y2;y3;y4;y5;y6;y7;y8];
FrequencyArray = zeros(FreqNum,1);
MagnitudeArray = zeros(FreqNum,1);
PhaseArray = zeros(FreqNum,1);
for i = 1:1:FreqNum
In_Signal = sys_identify_in(AllocationResult(i,2), AllocationResult(i,3):1:AllocationResult(i,4));
Out_Signal = sys_identify_out(AllocationResult(i,2), AllocationResult(i,3):1:AllocationResult(i,4));
InFFT = fft(In_Signal);
OutFFT = fft(Out_Signal);
FrequencyArray(i) = AllocationResult(i,1);
MagnitudeArray(i) = abs(OutFFT(2)/InFFT(2)); % 第一个信号为直流分量 第二个信号为主分量
PhaseArray(i) = angle(OutFFT(2)/InFFT(2))/pi*180; % 角度处理为角度制
end
obj = [FrequencyArray, MagnitudeArray, PhaseArray];
save(['MatlabData/','obj','.mat'],'obj');
ObjNum = 2;
ObjPath = {'MatlabData/obj.mat','MatlabData/obj1.mat'};
ObjLegend = {'test','test1'};
ObjTitle = '频域特性';
comparison(ObjNum,ObjPath,ObjLegend,ObjTitle);
运行结果如下：
可以翻一下，前面我们DSP代码中的测试代码就是输出是输入的10倍（20dB）相位相差pi/3(60°)。第一次测试有问题的原因数据解析的问题。
到此为止，DSP无通信接口系统辨识就结束了。 相关代码我会上传备份一手。
DSP无通信接口高效频域辨识-嵌入式文档类资源-CSDN文库
下面给一个直观的辨识视频，舒服
高速电机FOC iq→iq 电流环被控对象辨识
感谢您的阅读，欢迎留言讨论、收藏、点赞、分享。
扫一扫
3786
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
