]> Creatis software - gdcm.git/blob - src/gdcmmpeg2/src/mpeg2enc/fdctref.c
COMP: A few stupid cast needed for vs7
[gdcm.git] / src / gdcmmpeg2 / src / mpeg2enc / fdctref.c
1 /* fdctref.c, forward discrete cosine transform, double precision           */
2
3 /* Copyright (C) 1996, MPEG Software Simulation Group. All Rights Reserved. */
4
5 /*
6  * Disclaimer of Warranty
7  *
8  * These software programs are available to the user without any license fee or
9  * royalty on an "as is" basis.  The MPEG Software Simulation Group disclaims
10  * any and all warranties, whether express, implied, or statuary, including any
11  * implied warranties or merchantability or of fitness for a particular
12  * purpose.  In no event shall the copyright-holder be liable for any
13  * incidental, punitive, or consequential damages of any kind whatsoever
14  * arising from the use of these programs.
15  *
16  * This disclaimer of warranty extends to the user of these programs and user's
17  * customers, employees, agents, transferees, successors, and assigns.
18  *
19  * The MPEG Software Simulation Group does not represent or warrant that the
20  * programs furnished hereunder are free of infringement of any third-party
21  * patents.
22  *
23  * Commercial implementations of MPEG-1 and MPEG-2 video, including shareware,
24  * are subject to royalty fees to patent holders.  Many of these patents are
25  * general enough such that they are unavoidable regardless of implementation
26  * design.
27  *
28  */
29
30 #include <math.h>
31
32 #include "config.h"
33
34 #ifndef PI
35 # ifdef M_PI
36 #  define PI M_PI
37 # else
38 #  define PI 3.14159265358979323846
39 # endif
40 #endif
41
42 /* global declarations */
43 void init_fdct _ANSI_ARGS_((void));
44 void fdct _ANSI_ARGS_((short *block));
45
46 /* private data */
47 static double c[8][8]; /* transform coefficients */
48
49 void init_fdct()
50 {
51   int i, j;
52   double s;
53
54   for (i=0; i<8; i++)
55   {
56     s = (i==0) ? sqrt(0.125) : 0.5;
57
58     for (j=0; j<8; j++)
59       c[i][j] = s * cos((PI/8.0)*i*(j+0.5));
60   }
61 }
62
63 void fdct(block)
64 short *block;
65 {
66   int i, j, k;
67   double s;
68   double tmp[64];
69
70   for (i=0; i<8; i++)
71     for (j=0; j<8; j++)
72     {
73       s = 0.0;
74
75       for (k=0; k<8; k++)
76         s += c[j][k] * block[8*i+k];
77
78       tmp[8*i+j] = s;
79     }
80
81   for (j=0; j<8; j++)
82     for (i=0; i<8; i++)
83     {
84       s = 0.0;
85
86       for (k=0; k<8; k++)
87         s += c[i][k] * tmp[8*k+j];
88
89       block[8*i+j] = (int)floor(s+0.499999);
90       /*
91        * reason for adding 0.499999 instead of 0.5:
92        * s is quite often x.5 (at least for i and/or j = 0 or 4)
93        * and setting the rounding threshold exactly to 0.5 leads to an
94        * extremely high arithmetic implementation dependency of the result;
95        * s being between x.5 and x.500001 (which is now incorrectly rounded
96        * downwards instead of upwards) is assumed to occur less often
97        * (if at all)
98        */
99     }
100 }