DCTNew

        public static void Main(string[] args)
        {
            InitializeDCT();
            double[] ValIn = new double[] {
                -76, -73, -67, -62, -58, -67, -64, -55,
                -65, -69, -73, -38, -19, -43, -59, -56,
                -66, -69, -60, -15,  16, -24, -62, -55,
                -65, -70, -57, - 6,  26, -22, -58, -59,
                -61, -67, -60, -24, - 2, -40, -60, -58,
                -49, -63, -68, -58, -51, -60, -70, -53,
                -43, -57, -64, -69, -73, -67, -63, -45,
                -41, -49, -59, -60, -63, -52, -50, -34
            };
            double[] ValDCT = DCT(ValIn);
            double[] ValOut = IDCT(ValDCT);
            double[] ValDiff = new double[64];
            for (int i = 0; i < 64; i++) ValDiff[i] = Math.Round(ValIn[i] - ValOut[i]);
        }

        unsafe static double D( ulong val) => *(double*)&val;
        unsafe static  ulong I(double val) => *( ulong*)&val;

        static double[] Coeff;
        static double[] C;

        static void InitializeDCT(byte N = 8)
        {
            byte i, j;
            Coeff = new double[N * N];
            C     = new double[N    ];

            for (i = 0; i < N; i++)
            {
                if (i == 0) C[i] = 1 / Math.Sqrt(2);
                else        C[i] = 1;
                for (j = 0; j < N; j++)
                    Coeff[i * N + j] = 0.5 * C[i] * Math.Cos((2 * j + 1) * i / (N * 2.0) * Math.PI);
            }
        }

        static double[] DCT(double[] In, byte N = 8)
        {
            byte u, v, x, y;
            double[] Out = new double[N * N];

            for (v = 0; v < N; v++)
                for (u = 0; u < N; u++)
                    for (y = 0, Out[v * N + u] = 0; y < N; y++)
                        for (x = 0; x < N; x++)
                            Out[v * N + u] += In[y * N + x] * Coeff[v * N + y] * Coeff[u * N + x];
            return Out;
        }

        static double[] IDCT(double[] In, byte N = 8)
        {
            byte u, v, x, y;
            double[] Out = new double[N * N];

            for (v = 0; v < N; v++)
                for (u = 0; u < N; u++)
                    for (y = 0, Out[v * 3 + u] = 0; y < N; y++)
                        for (x = 0; x < N; x++)
                            Out[v * N + u] += In[y * N + x] * Coeff[y * N + v] * Coeff[x * N + u];
            return Out;
        }