C#计算功率谱密度

在C#中计算功率谱密度，可以使用FFT（快速傅里叶变换）来实现。下面是一个简单的示例代码：

using System;
using System.Numerics;

namespace PowerSpectralDensity
{
    class Program
    {
        static void Main(string[] args)
        {
            // 读取音频数据
            double[] data = ReadAudioData("audio.wav");

            // 计算FFT
            Complex[] fft = new Complex[data.Length];
            for (int i = 0; i < data.Length; i++)
            {
                fft[i] = new Complex(data[i], 0);
            }
            FourierTransform.FFT(fft, FourierTransform.Direction.Forward);

            // 计算功率谱密度
            double[] powerSpectrum = new double[data.Length / 2];
            for (int i = 0; i < data.Length / 2; i++)
            {
                powerSpectrum[i] = Math.Pow(fft[i].Magnitude, 2) / (data.Length / 2);
            }

            // 输出结果
            for (int i = 0; i < powerSpectrum.Length; i++)
            {
                Console.WriteLine("{0} Hz: {1} dB", i * (44100.0 / data.Length), 10 * Math.Log10(powerSpectrum[i]));
            }
        }

        static double[] ReadAudioData(string filename)
        {
            // TODO: 从音频文件中读取数据
            return new double[0];
        }
    }

    static class FourierTransform
    {
        public enum Direction
        {
            Forward,
            Backward
        }

        public static void FFT(Complex[] data, Direction direction)
        {
            int n = data.Length;
            int m = (int)Math.Log(n, 2);

            // 交换数据
            for (int i = 0; i < n; i++)
            {
                int j = ReverseBits(i, m);
                if (j > i)
                {
                    Complex temp = data[i];
                    data[i] = data[j];
                    data[j] = temp;
                }
            }

            // 计算FFT
            for (int s = 1; s <= m; s++)
            {
                int m2 = 1 << s;
                Complex w = new Complex(1, 0);
                Complex wm = new Complex(Math.Cos(Math.PI / (m2 / 2)), Math.Sin(Math.PI / (m2 / 2)));
                if (direction == Direction.Backward)
                {
                    wm = Complex.Conjugate(wm);
                }

                for (int j = 0; j < m2 / 2; j++)
                {
                    for (int k = j; k < n; k += m2)
                    {
                        Complex t = w * data[k + m2 / 2];
                        Complex u = data[k];
                        data[k] = u + t;
                        data[k + m2 / 2] = u - t;
                    }

                    w *= wm;
                }
            }

            if (direction == Direction.Forward)
            {
                for (int i = 0; i < n; i++)
                {
                    data[i] /= n;
                }
            }
        }

        static int ReverseBits(int value, int bits)
        {
            int result = 0;
            for (int i = 0; i < bits; i++)
            {
                result <<= 1;
                result |= value & 1;
                value >>= 1;
            }
            return result;
        }
    }
}

在上面的代码中，我们首先读取音频数据，然后计算FFT，并通过FFT计算功率谱密度。最后输出结果，其中频率为x Hz的功率谱密度为y dB。需要注意的是，上面的代码仅提供了一个简单的示例，实际使用时需要根据具体情况进行修改。