使用 fopen_s读取一个wav文件并对其进行fft和STFT并分别输出结果c语言代码不使用fftw库函数并解释代码

以下是使用 C 语言代码实现读取 WAV 文件并进行 FFT 和 STFT 的示例代码：

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

#define PI 3.14159265358979323846

// WAV 文件头结构体
typedef struct {
    char chunkID[4];
    unsigned int chunkSize;
    char format[4];
    char subchunk1ID[4];
    unsigned int subchunk1Size;
    unsigned short audioFormat;
    unsigned short numChannels;
    unsigned int sampleRate;
    unsigned int byteRate;
    unsigned short blockAlign;
    unsigned short bitsPerSample;
    char subchunk2ID[4];
    unsigned int subchunk2Size;
} WAV_HEADER;

// 计算幅度谱
void amplitudeSpectrum(int N, double *x, double *y) {
    int i;
    for (i = 0; i < N; ++i) {
        y[i] = sqrt(x[i] * x[i] + y[i] * y[i]);
    }
}

// 计算FFT
void fft(int N, double *x, double *y) {
    int i, j, k;
    double c, s, t1, t2;

    // 按位倒置
    j = 0;
    for (i = 0; i < N - 1; ++i) {
        if (i < j) {
            t1 = x[i];
            x[i] = x[j];
            x[j] = t1;
            t1 = y[i];
            y[i] = y[j];
            y[j] = t1;
        }
        k = N / 2;
        while (k <= j) {
            j -= k;
            k /= 2;
        }
        j += k;
    }

    // FFT 计算
    int m = log2(N);
    for (i = 0; i < m; ++i) {
        int n = pow(2, i);
        for (j = 0; j < N; j += 2 * n) {
            for (k = 0; k < n; ++k) {
                c = cos(-PI * k / n);
                s = sin(-PI * k / n);
                t1 = c * x[j + n + k] - s * y[j + n + k];
                t2 = s * x[j + n + k] + c * y[j + n + k];
                x[j + n + k] = x[j + k] - t1;
                y[j + n + k] = y[j + k] - t2;
                x[j + k] += t1;
                y[j + k] += t2;
            }
        }
    }
}

// 计算STFT
void stft(int N, int M, double *x, double *y, double *X, double *Y) {
    int i, j, k;
    double *win = (double *) malloc(sizeof(double) * M);
    double *tmpx = (double *) malloc(sizeof(double) * N);
    double *tmpy = (double *) malloc(sizeof(double) * N);

    // 窗函数
    for (i = 0; i < M; ++i) {
        win[i] = 0.5 - 0.5 * cos(2 * PI * i / (M - 1));
    }

    // STFT 计算
    for (i = 0; i < N - M; i += M) {
        for (j = 0; j < M; ++j) {
            tmpx[j] = x[i + j] * win[j];
            tmpy[j] = y[i + j] * win[j];
        }
        fft(M, tmpx, tmpy);
        for (j = 0; j < M / 2; ++j) {
            X[i / M * (M / 2) + j] = tmpx[j];
            Y[i / M * (M / 2) + j] = tmpy[j];
        }
    }

    free(win);
    free(tmpx);
    free(tmpy);
}

int main() {
    FILE *fp;
    WAV_HEADER header;
    double *data;
    double *x;
    double *y;
    double *X;
    double *Y;
    int N, M;
    int i;

    // 读取 WAV 文件头
    fp = fopen("test.wav", "rb");
    if (!fp) {
        printf("Failed to open file.\n");
        exit(1);
    }
    fread(&header, sizeof(WAV_HEADER), 1, fp);

    // 计算数据长度
    N = header.subchunk2Size / (header.bitsPerSample / 8) / header.numChannels;
    M = 1024;  // STFT 窗长

    // 分配内存
    data = (double *) malloc(sizeof(double) * N);
    x = (double *) malloc(sizeof(double) * N);
    y = (double *) malloc(sizeof(double) * N);
    X = (double *) malloc(sizeof(double) * N / M * M / 2);
    Y = (double *) malloc(sizeof(double) * N / M * M / 2);

    // 读取数据
    for (i = 0; i < N; ++i) {
        data[i] = 0;
        fread(&data[i], header.bitsPerSample / 8, 1, fp);
    }

    // 分离通道
    for (i = 0; i < N; ++i) {
        if (header.numChannels == 1) {
            x[i] = data[i];
            y[i] = 0;
        } else {
            x[i] = data[i * 2];
            y[i] = data[i * 2 + 1];
        }
    }

    // 计算 STFT
    stft(N, M, x, y, X, Y);

    // 输出结果
    for (i = 0; i < N / M * M / 2; ++i) {
        printf("%lf\n", X[i]);
    }

    // 释放内存
    free(data);
    free(x);
    free(y);
    free(X);
    free(Y);

    return 0;
}

在上面的代码中，我们定义了一个 WAV 文件头结构体 `WAV_HEADER`，用于读取 WAV 文件头信息。然后我们定义了三个函数：`amplitudeSpectrum` 用于计算幅度谱，`fft` 用于计算 FFT，`stft` 用于计算 STFT。其中，FFT 的实现采用了蝴蝶算法，STFT 的实现采用了重叠窗口法。最后，在 `main` 函数中，我们读取 WAV 文件数据，将数据分离通道，并计算 STFT，并输出结果。

需要注意的是，本代码中的 FFT 和 STFT 实现仅用于教学和学习，实际应用中建议使用专业的 FFT 和 STFT 库函数，如 FFTW 库函数。