C语言信号处理，离线版（全域、后处理）中值滤波和低通滤波
一些说明&简要介绍
对于一些数据，需要在pc或者云端进行滤波处理，比如这里的中值（去基线）和低通（去毛刺）。因为平台强大，当数据量大时，如果滤波做的不是太复杂，或者说做的比较简单，那么处理起来还是比较快的，那么就能实用。
这里有个例子，因为是后处理，则信号处理前是多长（多少点），处理后就是多少点。对于中值滤波，需要将边沿点进行扩展（补点），或者在边沿进行窗口缩小处理。对于低通滤波，则确定需要对边沿进行补点，并合理处理延时。
简单来说，中值滤波边沿补点处理有很多方式，matlab的有中值滤波系统函数median，但是不带这个参数（印象中好像是带的，今天再次去查了很久，发现没有，不知道是不是版本问题还是其他），直接给答案吧：常见的可以是补零、边沿值复制填充和对称填充。这个对称填充，又可以是一次对称填充和逐元素对称填充，等等。肯定还有很多其他的填充方式，暂且不深究。我们这里采用补零填充。
低通滤波，我们这里用FIR，所以是固定延时的（延时点数：阶数/2），延时可以解释下：即滤到当前点时，结果是当前点前阶数/2点的结果。所以，要滤出当前点的值，则需要滤到当前点后阶数/2点，其实，就是将滤波器系数窗口（FIR滤波器系数窗口是一个对称的类似汉宁窗的窗口）的中间值（我们选择偶数阶的滤波器，那么系数中间值是一个值，其两边的系数刚好对称）对齐需要滤波的信号点，滤出的就是这个点无延时的结果。如此说来，也就知道了我们的滤波是怎么处理的，即对称扩展。需要注意的是，不是一次扩展，而是对边沿的每个点都对称扩展，即：因为其前面点不够，将后面的阶数/2个点对称扩展到该点的前面，然后用滤波窗口滤波，实现了无延时且滤波正确的结果。
比较有意思的是，我们既可以单独测试中值和低通滤波，还可以将两个滤波器串联起来——其实就是两次滤波，因为是后处理，信号可以处理成无延时，也就是说，只要滤波正确，多少次滤波都可以。
下面是C代码：
main.c
//#include<stdio.h>
//#include<stdlib.h>
//#include<math.h>
//
//
//#include"Median.h"
//
//
//
//#define M_PI 3.1415
//#define N 200
//#define F1 3
//#define F2 40
//
//int main()
//{
// // 生成正弦信号
// int a[N] = { 0 }, b[N] = { 0 };
// int magnify = 100;
// for (int i = 0; i < N; i++)
// {
// a[i] = (int)(magnify *sin(2 * M_PI * i * F1 / N) + 0.25 * magnify *sin(2 * M_PI * i * F2 / N));
// }
//
// FILE* fp;
//
// fp = fopen("text.txt", "w");
// if (feof(fp))
// {
// printf("NULL");
// exit(0); // 表示如果读取为空文件就正常退出
// }
//
// for (int i = 0; i < N; i++)
// fprintf(fp, "%d\n", a[i]);
// fclose(fp);
//
// int fir_len = 10;
// int m;
//
// // 中值滤波
// m = Median(a, b, N, fir_len);
//
//
// fp = fopen("text_fir.txt", "w");
// if (feof(fp))
// {
// printf("NULL");
// exit(0);//表示如果读取为空文件就正常退出
// }
//
// for (int i = 0; i < N; i++)
// fprintf(fp, "%d\n", b[i]);
// fclose(fp);
//
//
// return 0;
//}
//
#include<stdio.h>
#include<stdlib.h>
#include <string.h>
#include "Lowpass.h"
#include "Median.h"
#define SampleRate (256)
int main()
{
//errno_t err;
FILE* file;
//#pragma warning(disable:4996)
file = fopen("mydata.txt", "rb"); // 以二进制模式打开文件供读取
//err = fopen_s(&file, "D:\WorkDoc\algoPro\AlgoSDKPro\AlgoSDK\mydata.txt", "rb");
if (file == NULL)
{
printf("找不到文件");
return 0;
}
fseek(file, 0, SEEK_END); // 将文件指针移动到文件末尾
long fileSize = ftell(file); // 获取文件大小
printf("文件大小%ld", fileSize);
rewind(file); // 将文件指针重新定位到文件开头
char* buffer = (char*)malloc(fileSize + 1); // 分配足够大的缓冲区来存储文件内容
if (buffer == NULL) {
perror("Memory allocation failed");
return 1;
}
fread(buffer, 1, fileSize, file); // 一次性读取文件全部内容
buffer[fileSize] = '\0'; // 添加字符串结束符
//char buffer[100]; // 假设每次读取最多包含 100 个字符
int comma_count = 0; // 逗号计数器
//while (fgets(buffer, sizeof(buffer), file)) {
for (int i = 0; buffer[i] != '\0'; i++) {
if (buffer[i] == ',') {
comma_count++;
}
}
//}
printf("Number of commas in the file: %d\n", comma_count);
fclose(file);
if (comma_count > 0)
{
int data_size = comma_count + 1;
// 使用 malloc 动态分配数组内存空间
//int* data_in = (int*)malloc(data_size * sizeof(int));
//int* data_out = (int*)malloc(data_size * sizeof(int));
// 使用 calloc 动态分配数组内存空间并初始化为0
int* data_in = (int*)calloc(data_size, sizeof(int));
int* data_out = (int*)calloc(data_size, sizeof(int));
int* data_out_finish = (int*)calloc(data_size, sizeof(int));
//
//int data_in[data_size] = { 0 };
//int data_out[data_size] = { 0 };
char* token = strtok(buffer, ","); //strtok引用#include <string.h>
int index = 0;
while (token != NULL) {
//double number = atof(token); // 将分割出的字符串转换为浮点数
int number = atoi(token);
//printf("%s\n", token);
data_in[index++] = number;
token = strtok(NULL, ",");
}
printf("输出原始的ADC数值：");
//输出滤波前的数值
for (int j = 0; j < data_size; j++)
{
printf("%d ", data_in[j]);
}
调用算法 start
LowpassInit();
int lRet = Lowpass(data_in, data_size, data_out, SampleRate);
if (lRet) {
free(buffer);
return 1;
}
int fir_len = 10;
int mRet = Median(data_out, data_out_finish, data_size, fir_len);
调用算法 end
printf("\n------------------------------------------------------------------\n");
printf(" \n输出滤波后的数值");
//输出滤波后的数值
for (int j = 0; j < data_size; j++)
{
printf("%d ", data_out[j]);
}
}
free(buffer);
return 0;
}
开头的注释部分，是C自己生成正弦型信号测试中值滤波的。后面部分是读入数据文件测试低通和中值串联的。
Lowpass.h
#pragma once
extern void LowpassInit();
extern int Lowpass(int* data_in,int data_len, int* data_out,int SampleRate);Lowpass.c
#include"stdio.h"
#include"stdlib.h"
#include"string.h"
#include "assert.h"
// FIR 45Hz低通
static float fir_coef_b[129] ={
0.001364, -0.001009, -0.001358, -0.000845, 0.000328, 0.000870, 0.000046, -0.001164, -0.001094, 0.000425, 0.001572, 0.000723, -0.001273, -0.001857, -0.000006, 0.002128,
0.001681, -0.001116, -0.002802, -0.000929, 0.002459, 0.002985, -0.000451, -0.003734, -0.002413, 0.002331, 0.004545, 0.000925, -0.004396, -0.004473, 0.001453, 0.006163,
0.003177, -0.004447, -0.007040, -0.000497, 0.007527, 0.006437, -0.003443, -0.009946, -0.003896, 0.008190, 0.010827, -0.000795, -0.012942, -0.009279, 0.007504, 0.016601,
0.004484, -0.015728, -0.017795, 0.004292, 0.024643, 0.014802, -0.017997, -0.033204, -0.004871, 0.039100, 0.040318, -0.019479, -0.078933, -0.045026, 0.108777, 0.295845,
0.380006, 0.295845, 0.108777, -0.045026, -0.078933, -0.019479, 0.040318, 0.039100, -0.004871, -0.033204, -0.017997, 0.014802, 0.024643, 0.004292, -0.017795, -0.015728,
0.004484, 0.016601, 0.007504, -0.009279, -0.012942, -0.000795, 0.010827, 0.008190, -0.003896, -0.009946, -0.003443, 0.006437, 0.007527, -0.000497, -0.007040, -0.004447,
0.003177, 0.006163, 0.001453, -0.004473, -0.004396, 0.000925, 0.004545, 0.002331, -0.002413, -0.003734, -0.000451, 0.002985, 0.002459, -0.000929, -0.002802, -0.001116,
0.001681, 0.002128, -0.000006, -0.001857, -0.001273, 0.000723, 0.001572, 0.000425, -0.001094, -0.001164, 0.000046, 0.000870, 0.000328, -0.000845, -0.001358, -0.001009, 0.001364,
};
static float gain = 1.01;
typedef struct {
int DataBuf[129];
float fDataBuf[129];
float DataAfFir;
}LowpassPra_t;
LowpassPra_t LowpassPra = { 0 };
void fir_filter_zhh(float* sig_in, float* sig_out,int sig_len)
{
assert(sig_len > 0);
int i, j;
for (i = 0; i < sig_len; i++) {
sig_out[i] = 0;
for (j = 0; j < 129; j++) {
sig_out[i] += sig_in[i - j] * fir_coef_b[j];
}
sig_out[i] *= gain;
}
}
// 数组逆序 float版
void reverse_arry_float(float* arry,int arry_len) // arry_len可以为1，不报错，此时arry只有一个元素，调用后不变
{
int i = 0; //循环变量1, i的值为数组第一个元素的下标
int j = arry_len - 1; //循环变量2, j的值为数组最后一个元素的下标
float fir_idx_buf; //互换时的中间存储变量
for (; i < j; ++i, --j) /*因为i和j已经初始化过了, 所以表达式1可以省略, 但表达式1后面的分号不能省。*/
{
fir_idx_buf = arry[i];
arry[i] = arry[j];
arry[j] = fir_idx_buf;
}
}
void LowpassInit()
{
memset(LowpassPra.DataBuf, 0, sizeof(LowpassPra.DataBuf));
}
// 数据按帧传入，同时要传入帧号
int Lowpass(int* data_in,int data_len, int* data_out,int SampleRate)
{
// 限制采样率必须是256
if (SampleRate != 256) {
//printf("SampleRate error!\n");
return 1;
}
int i,j;
for (i = 0;i < data_len;i++) {
if( i < 64 ) { // 向左边扩展
memcpy(&LowpassPra.DataBuf[0], &data_in[i + 1], 64 * sizeof(data_in[0]));
memcpy( &LowpassPra.DataBuf[64], &data_in[i], 65 * sizeof(data_in[0]) );
reverse_arry_float(&LowpassPra.DataBuf[0], 64);
for (j = 0; j <= 128;j++) {
LowpassPra.fDataBuf[j] = (float)LowpassPra.DataBuf[j];
}
}
else if( i < data_len - 64 ) {
for (j = 0; j <= 128; j++) {
LowpassPra.fDataBuf[j] = (float)data_in[i - 64 + j];
}
}
else { // 向右边扩展
memcpy(&LowpassPra.DataBuf[0], &data_in[i - 64], 65 * sizeof(data_in[0]));
memcpy(&LowpassPra.DataBuf[64 + 1], &data_in[i - 64], 64 * sizeof(data_in[0]));
reverse_arry_float(&LowpassPra.DataBuf[64 + 1], 64);
for (j = 0; j <= 128; j++) {
LowpassPra.fDataBuf[j] = (float)LowpassPra.DataBuf[j];
}
}
fir_filter_zhh(&LowpassPra.fDataBuf[128], &LowpassPra.DataAfFir, 1);
data_out[i] = (int)LowpassPra.DataAfFir;
}
return 0;
}
Median.h
#pragma once
extern int Median(int* sig, int* sig_out,int sig_len,int fir_len);Median.c
#include<stdio.h>
#include<stdlib.h>
#include<string.h>
#include <assert.h>
void quick_sort_int(int* a, int low, int high);
中值滤波函数
// C代码方式模拟MATLAB中值滤波函数
// sig: 输入信号，必须是行向量，且元素个数至少是2
// sig_len：信号总长度
// fir_len : 窗口长度, 2到信号长度之间，或者后续再另外设置最小最大值。
int Median(int* sig, int* sig_out,int sig_len,int fir_len)
{
// 输入信号长度检查，必须大于1
assert(sig_len > 1);
// 检查fir_len，整数型，最少2，最大为信号sig的长度
assert(fir_len > 1 && fir_len <= sig_len);
int i;
int median_tmp;
// 需要缓存的数组，无论哪种扩展，无论fir_len奇偶，都是补fir_len - 1个点
int* sig_buf = (int*)malloc((sig_len + fir_len - 1) * sizeof(int));
if (sig_buf == NULL) {//判空
printf("sig_buf malloc error!\n");//打印错误信息
return 1;
}
memset(sig_buf, 0, (sig_len + fir_len - 1) * sizeof(int));
// 补零扩展
int* buf = (int*)malloc(fir_len * sizeof(int));
if (buf == NULL) {//判空
free(sig_buf);
printf("buf malloc 1 error!\n");//打印错误信息
return 1;
}
if (fir_len % 2 == 0) { // 偶数
memset(sig_buf, 0, (fir_len / 2) * sizeof(int));
memcpy(sig_buf + fir_len / 2, sig, sig_len * sizeof(int));
memset(sig_buf + fir_len / 2 + sig_len, 0, (fir_len / 2 - 1) * sizeof(int));
// 滤波过程
for (i = fir_len / 2; i < sig_len + fir_len / 2; i++) {
memcpy(buf, &sig_buf[i - fir_len / 2], fir_len * sizeof(int));
quick_sort_int(buf, 0, fir_len - 1);// 升序排序
median_tmp = (int)(buf[fir_len / 2 - 1] + buf[fir_len / 2]) / 2;
sig_out[i - fir_len / 2] = median_tmp;
}
}
else {// 奇数
memset(sig_buf, 0, ((fir_len - 1) / 2) * sizeof(int));
memcpy(sig_buf + (fir_len - 1) / 2, sig, sig_len * sizeof(int));
memset(sig_buf + (fir_len - 1) / 2 + sig_len, 0, ((fir_len - 1) / 2) * sizeof(int));
// 滤波过程
for (i = (fir_len - 1) / 2; i < sig_len + (fir_len - 1) / 2; i++) {
memcpy(buf, &sig_buf[i - (fir_len - 1) / 2], fir_len * sizeof(int));
quick_sort_int(buf, 0, fir_len - 1);// 升序排序
median_tmp = buf[(fir_len - 1) / 2];
sig_out[i - (fir_len - 1) / 2] = median_tmp;
}
}
free(sig_buf);
free(buf);
return 0;// 正常返回0
}
void quick_sort_int(int* a, int low, int high)
{
int i = low;//第一位
int j = high;//最后一位
int key = a[i]; //将第一个数作为基准值-- 先找到一个基准值
while (i < j)
{
while (i < j && a[j] >= key)
{
j--;
}
a[i] = a[j];
while (i < j && a[i] <= key)
{
i++;
}
a[j] = a[i];
}
a[i] = key;
if (i - 1 > low)
{
quick_sort_int(a, low, i - 1);
}
if (i + 1 < high)
{
quick_sort_int(a, i + 1, high);
}
}
原文地址：https://blog.csdn.net/weixin_40194697/article/details/142363622
免责声明：本站文章内容转载自网络资源，如本站内容侵犯了原著者的合法权益，可联系本站删除。更多内容请关注自学内容网（zxcms.com）！