Использование КИХ-фильтра
15.10.2015, 13:21. Показов 2865. Ответов 1
Здравствуйте. В матлабе создал сигнал, частота дискретизации 192 кГц, частоты сигнала 50 кГц и 75 кГц. Получил через fdatool коэффициенты (их 487) для полосового КИХ-фильтра (оставляет только 75 кГц). Сигнал и коэффициенты записал в соответствующие файлы. Попытался на С++ реализовать фильтрацию и посмотреть результат снова в матлабе, но видимо, что-то кардинально делаю не так, потому что результат совсем не тот, что нужен. Можете подсказать, что не так?
| Matlab M | 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
| Fs=192000; % частота дискретизации
t=1/Fs:1/Fs:1; % дискретное время
f1=50000; % первая частота модуляции
f2=75000; % вторая частота модуляции
s_M=cos(2*pi*f1*t)+cos(2*pi*f2*t); % модулирующий сигнал
s_M=s_M.*10000; %чтобы были целые значения, а не double
Num=Num.*1000000; %это из fdatool
fid = fopen('D:\\Sasha\\dd\\Num', 'wb'); % открытие файла на запись
if fid == -1 % проверка корректности открытия
error('File is not opened');
end
fwrite(fid, Num, 'int'); % запись матрицы в файл
fclose(fid); % закрытие файла
fid = fopen('D:\\Sasha\\dd\\Signal', 'wb'); % открытие файла на запись
if fid == -1 % проверка корректности открытия
error('File is not opened');
end
fwrite(fid, s_M, 'int'); % запись матрицы в файл
fclose(fid); % закрытие файла |
|
| C++ | 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
| #include <iostream>
#include <math.h>
#include <fstream>
#include <stdio.h>
using namespace std;
int main()
{
int nFileLen;
int mFileLen;
int buffer[192000];
int result[192000];
int h[487];
FILE *f;
f = fopen("D:\\Sasha\\dd\\Signal", "rb"); //open file
if (f == NULL)
cout << "ошибка"; //error
else { //all right
fseek(f, 0, SEEK_END); //pointer to the end of file
nFileLen = ftell(f); //number of last byte
fseek(f, 0, SEEK_SET); //pointer to the start of file
for (int j = 0; j < nFileLen / (sizeof(int)); j++) { //read queue items
fread(&buffer[j], sizeof(int), 1, f);
};
fclose(f); //close file
};
FILE *p;
p = fopen("D:\\Sasha\\dd\\Num", "rb"); //open file
if (p == NULL)
cout << "ошибка"; //error
else { //all right
fseek(p, 0, SEEK_END); //pointer to the end of file
mFileLen = ftell(p); //number of last byte
fseek(p, 0, SEEK_SET); //pointer to the start of file
for (int j = 0; j < mFileLen / (sizeof(int)); j++) { //read queue items
fread(&h[j], sizeof(int), 1, p);
};
fclose(p); //close file
};
for (int i=0; i<192000; i++)
{
result[i]=0.;
for (int j=0; j<486; j++)// формула фильтра
{
result[i]+= h[j]*buffer[i-j];
}
}
FILE *lol;
lol = fopen("D:\\Sasha\\dd\\result", "wb");
for (int q = 0; q < 192000; q++) {
fwrite(&result[q], sizeof(int), 1, lol);
};
fclose(lol);
} |
|
| Matlab M | 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
| fid = fopen('D:\\Sasha\\dd\\result', 'rb'); % открытие файла на чтение
if fid == -1 % проверка корректности открытия
error('File is not opened');
end
B = fread(fid, 192000, 'int'); % чтение значений
B=B';
fclose(fid);% закрытие файла
FftL=Fs;
FftS=abs(fft(B,FftL));% Амплитуды преобразования Фурье сигнала
FftS=2*FftS./FftL;% Нормировка спектра по амплитуде
F=0:1:FftL/2-1/FftL;% Массив частот вычисляемого спектра Фурье
figure;
plot(F,FftS(1:length(F)));% Построение спектра Фурье сигнала
grid on; |
|
На форуме нашел реализацию матлабовской функции filter на c++ от vital792. Попробовал её применить, изменил только, что входной и выходной массив будут типа int, а не double и т.к. использую КИХ-фильтр, то вместо коэффициентов знаменателя передавал 1 0 0 0 ... Результат опять почему-то не тот 
| C++ | 1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
| #include <iostream>
#include <fstream>
#include <stdio.h>
#include <cmath>
#include <cstring>
using namespace std;
void filter(int ord, const double a[], const double b[], int length, int *x, int *y)
{
double *t = new double[ord+1];
memset(t, 0, (ord+1)*sizeof(double));
for(int i=0; i<length; i++)
{
double xi = x[i];
t[0] = b[0]*xi + t[1];
double t0 = t[0];
for(int j=1; j<ord; j++)
{
t[j] = b[j]*xi - a[j]*t0 + t[j+1];
}
t[ord] = b[ord]*xi - a[ord]*t0;
y[i] = a[0]*t0;
}
delete [] t;
}
int main()
{
int nFileLen;
int mFileLen;
int buffer[192000];
int result[192000];
double h[487];
double d[487];
FILE *f;
f = fopen("D:\\Sasha\\dd\\Signal", "rb"); //open file
if (f == NULL)
cout << "ошибка"; //error
else { //all right
fseek(f, 0, SEEK_END); //pointer to the end of file
nFileLen = ftell(f); //number of last byte
fseek(f, 0, SEEK_SET); //pointer to the start of file
for (int j = 0; j < nFileLen / (sizeof(int)); j++) { //read queue items
fread(&buffer[j], sizeof(int), 1, f);
};
fclose(f); //close file
};
FILE *p;
p = fopen("D:\\Sasha\\dd\\Num", "rb"); //open file
if (p == NULL)
cout << "ошибка"; //error
else { //all right
fseek(p, 0, SEEK_END); //pointer to the end of file
mFileLen = ftell(p); //number of last byte
fseek(p, 0, SEEK_SET); //pointer to the start of file
for (int j = 0; j < mFileLen / (sizeof(double)); j++) { //read queue items
fread(&h[j], sizeof(double), 1, p);
};
fclose(p); //close file
};
int order = 486;
memset(d, 0, (order+1)*sizeof(double));
d[0]=1;
filter(order, h, d, 192000, buffer, result);
FILE *lol;
lol = fopen("D:\\Sasha\\dd\\result", "wb");
for (int q = 0; q < 192000; q++) {
fwrite(&result[q], sizeof(int), 1, lol);
};
fclose(lol);
} |
|
Как я понял, ошибка была в задании коэффициентов фильтра
Добавлено через 1 минуту
Как правильно получить коэффициенты числителя и знаменателя передаточной функции из fdatool?
0
|
(
