FIR фильтр нижних частот - PullRequest
0 голосов
/ 30 марта 2020

Я пытаюсь реализовать низкочастотный фильтр FIR с окнами, но результат не соответствует.

Вот мой код:

#define N_FIR 51
main() {
    FILE* fin;  
    FILE* fout = NULL;  
    short *data16bits = NULL;

    fin = fopen("file3.wav", "rb");
    fread(&wavHeader, sizeof(wav_hdr), 1, fin);

    data16bits = (short*)malloc(wavHeader.subChunk2Size * sizeof(short));   
    fread(data16bits, wavHeader.bitsPerSample, wavHeader.subChunk2Size / 2, fin); 
    fclose(fin);   

    fout = fopen("FIRfilterLowPass.wav", "wb");
    fwrite(&wavHeader, sizeof(wav_hdr), 1, fout);

    double *fltr = (double *)malloc(N_FIR * sizeof(double));
    FIRfilterLowPass(1000, wavHeader.sampleRate, fltr);

    for (size_t i = 0; i < wavHeader.subChunk2Size; i++) {
        double result = 0;   
        short result_short = 0;
        for (size_t j = 0; j < N_FIR && (i - j) >= 0; j++) {
            result += fltr[j] * data16bits[i - j]; 
        }
        result_short = (short)result;   
        fwrite(&result_short, sizeof(short), 1, fout);
    }
}

void FIRfilterLowPass(double f_c, double f_samp, double* fltr) {
    f_c = f_c / f_samp;
    double w_c = 2 * PI * f_c;
    int middle = N_FIR / 2;

    for (int n = -N_FIR / 2; n <= N_FIR / 2; n++) {
        if (n == 0) {
            fltr[middle] = 2 * f_c;
        }
        else {
            fltr[middle + n] = sin(2 * PI * f_c * n) / (PI * n);
        }
    }

    for (int n = 0; n < N_FIR; n++) {
        fltr[n] = fltr[n] * (0.5f + 0.5f * cos((2 * PI * n) / N_FIR));
    }
}

Что я делаю не так, потому что волна выглядит так:

enter image description here

Но она должна выглядеть примерно так:

enter image description here

Спасибо за ваше время!

...