音频降风躁算法C代码

以下是一种基于数字信号处理的音频降风躁算法的C代码：

#include <stdio.h> #include <stdlib.h> #include <math.h>

#define PI 3.14159265358979323846

// 定义滤波器系数 float h[61] = {0.00000000, -0.00002213, -0.00006332, -0.00011826, -0.00018150, -0.00024461, -0.00029729, -0.00032759, -0.00032304, -0.00027194, -0.00016231, 0.00000000, 0.00019219, 0.00037248, 0.00050701, 0.00056778, 0.00054344, 0.00044728, 0.00031479, 0.00018493, 0.00008977, 0.00004031, 0.00000739, -0.00000000, -0.00000739, -0.00004031, -0.00008977, -0.00018493, -0.00031479, -0.00044728, -0.00054344, -0.00056778, -0.00050701, -0.00037248, -0.00019219, -0.00000000, 0.00016231, 0.00027194, 0.00032304, 0.00032759, 0.00029729, 0.00024461, 0.00018150, 0.00011826, 0.00006332, 0.00002213, 0.00000000, -0.00001135, -0.00001744, -0.00001714, -0.00001328, -0.00000824, -0.00000343, 0.00000000, 0.00000123, 0.00000136, 0.00000107, 0.00000069, 0.00000041};

// 定义音频帧大小和采样率 #define FRAME_SIZE 512 #define SAMPLE_RATE 16000

// 定义FFT所需的结构体 typedef struct { float real; float imag; } Complex;

Complex FFT(Complex x, Complex w) { Complex y; y.real = x.real * w.real - x.imag * w.imag; y.imag = x.imag * w.real + x.real * w.imag; return y; }

void FFT_Recursive(Complex* x, int N) { if (N <= 1) return;

Complex *even = (Complex*) malloc(N/2 * sizeof(Complex));
Complex *odd = (Complex*) malloc(N/2 * sizeof(Complex));

for (int i = 0; i < N/2; i++) {
    even[i] = x[2*i];
    odd[i] = x[2*i+1];
}

FFT_Recursive(even, N/2);
FFT_Recursive(odd, N/2);

for (int i = 0; i < N/2; i++) {
    Complex w;
    w.real = cos(2*PI*i/N);
    w.imag = -sin(2*PI*i/N);
    x[i] = FFT(even[i], w);
    x[i+N/2] = FFT(odd[i], w);
}

free(even);
free(odd);

}

void FFT_Iterative(Complex* x, int N) { int i, j, k, n1, n2, a; float c, s, t1, t2; Complex tx, t;

// Bit-reverse
j = 0;
n2 = N / 2;
for (i = 1; i < N - 1; i++) {
    n1 = n2;
    while (j >= n1) {
        j = j - n1;
        n1 = n1 / 2;
    }
    j = j + n1;

    if (i < j) {
        tx = x[i];
        x[i] = x[j];
        x[j] = tx;
    }
}

// FFT
n1 = 0;
n2 = 1;

for (i = 0; i < log2(N); i++) {
    n1 = n2;
    n2 = n2 + n2;
    a = 0;

    for (j = 0; j < n1; j++) {
        c = cos(2 * PI * a / n2);
        s = -sin(2 * PI * a / n2);
        a += 1 << (log2(N) - i - 1);

        for (k = j; k < N; k += n2) {
            t = FFT(x[k + n1], (Complex) {c, s});
            tx = x[k];
            x[k] = (Complex) {tx.real + t.real, tx.imag + t.imag};
            x[k + n1] = (Complex) {tx.real - t.real, tx.imag - t.imag};
        }
    }
}

}

void FFT_Compute(Complex* x, int N, int inverse) { if (inverse) { for (int i = 0; i < N; i++) { x[i].imag = -x[i].imag; } }

FFT_Iterative(x, N);

if (inverse) {
    for (int i = 0; i < N; i++) {
        x[i].real = x[i].real / N;
        x[i].imag = -x[i].imag / N;
    }
}

}

void Filter(Complex* x, int N) { int i; Complex y[FRAME_SIZE];

// 将实数数组转换为复数数组
for (i = 0; i < N; i++) {
    y[i].real = x[i];
    y[i].imag = 0.0;
}

// 对复数数组进行FFT
FFT_Compute(y, N, 0);

// 将滤波器系数也转换为复数数组
Complex h_c[FRAME_SIZE];
for (i = 0; i < 61; i++) {
    h_c[i].real = h[i];
    h_c[i].imag = 0.0;
}

// 对滤波器系数进行FFT
FFT_Compute(h_c, N, 0);

// 对频域信号进行滤波
for (i = 0; i < N; i++) {
    y[i] = FFT(y[i], h_c[i]);
}

// 对滤波后的频域信号进行IFFT
FFT_Compute(y, N, 1);

// 将复数数组转换为实数数组
for (i = 0; i < N; i++) {
    x[i] = y[i].real;
}

}

void Denoise(float* audio, int audio_len) { int i, j, nframes = audio_len / FRAME_SIZE;

for (i = 0; i < nframes; i++) {
    // 从音频中提取一帧
    float frame[FRAME_SIZE];
    for (j = 0; j < FRAME_SIZE; j++) {
        frame[j] = audio[i*FRAME_SIZE+j];
    }

    // 应用滤波器进行降噪
    Complex x[FRAME_SIZE];
    for (j = 0; j < FRAME_SIZE; j++) {
        x[j].real = frame[j];
        x[j].imag = 0.0;
    }
    Filter(x, FRAME_SIZE);
    for (j = 0; j < FRAME_SIZE; j++) {
        frame[j] = x[j].real;
    }

    // 将处理后的帧写回音频中
    for (j = 0; j < FRAME_SIZE; j++) {
        audio[i*FRAME_SIZE+j] = frame[j];
    }
}

}

int main() { // 读取音频文件 FILE* file = fopen("audio.wav", "rb"); if (file == NULL) { printf("Failed to open audio file."); return 1; }

// 解析音频文件头
char chunk_id[4], format[4], subchunk1_id[4], subchunk2_id[4];
int chunk_size, format_chunk_size, subchunk1_size, audio_format, num_channels, sample_rate, byte_rate, bits_per_sample, subchunk2_size;
fread(chunk_id, 1, 4, file);
fread(&chunk_size, 4, 1, file);
fread(format, 1, 4, file);
fread(subchunk1_id, 1, 4, file);
fread(&subchunk1_size, 4, 1, file);
fread(&audio_format, 2, 1, file);
fread(&num_channels, 2, 1, file);
fread(&sample_rate, 4, 1, file);
fread(&byte_rate, 4, 1, file);
fread(&bits_per_sample, 2, 1, file);
fread(subchunk2_id, 1, 4, file);
fread(&subchunk2_size, 4, 1, file);

// 读取音频数据
int audio_len = subchunk2_size / (bits_per_sample / 8);
float* audio = (float*) malloc(audio_len * sizeof(float));

if (bits_per_sample == 16) {
    short s;
    for (int i = 0; i < audio_len; i++) {
        fread(&s, 2, 1, file);
        audio[i] = (float) s / 32768.0;
    }
} else if (bits_per_sample == 32) {
    int i;
    for (i = 0; i < audio_len; i++) {
        fread(&audio[i], 4, 1, file);
    }
} else {
    printf("Unsupported audio format. Only 16-bit and 32-bit PCM are supported.");
    return 1;
}

// 关闭音频文件
fclose(file);

// 进行降风躁处理
Denoise(audio, audio_len);

// 将处理后的音频写回文件
file = fopen("audio_denoised.wav", "wb");
if (file == NULL) {
    printf("Failed to open output audio file.");
    return 1;
}

fwrite(chunk_id, 1, 4, file);
fwrite(&chunk_size, 4, 1, file);
fwrite(format, 1, 4, file);
fwrite(subchunk1_id, 1, 4, file);
fwrite(&subchunk1_size, 4, 1, file);
fwrite(&audio_format, 2, 1, file);
fwrite(&num_channels, 2, 1, file);
fwrite(&sample_rate, 4, 1, file);
fwrite(&byte_rate, 4, 1, file);
fwrite(&bits_per_sample, 2, 1, file);
fwrite(subchunk2_id, 1, 4, file);
fwrite(&subchunk2_size, 4, 1, file);

if (bits_per_sample == 16) {
    short s;
    for (int i = 0; i < audio_len; i++) {
        s = (short) (audio[i] * 32768.0);
        fwrite(&s, 2, 1, file);
    }
} else if (bits_per_sample == 32) {
    int i;
    for (i = 0; i < audio_len; i++) {
        fwrite(&audio[i], 4, 1, file);
    }
}

fclose(file);

return 0;