# Compute the Discrete Cosine Transform

Implement the Discrete Cosine Transform (DCT). This may implemented as either a function or a program and the sequence can be given as either an argument or using standard input. Your program must be able to handle sequences of any length (assuming a hypothetical version of your language which has no limits on things like memory and integer size).

There was a previous challenge for DFT, now lets compute DCT! ;-)

Use the DCT-II definition from Wikipedia: Your program takes a sequence xn as input, and must produce the corresponding sequence Xk. In this formula, the cosine is in radians.

## Rules

• This is so the shortest solution wins.
• Builtins that compute the DCT in forward or backward (also known as inverse) directions are not allowed.
• Floating-point inaccuracies will not be counted against you.
• You don't have to use the exact algorithm shown here, as long as you produce the same results.
• Mind including an explanation on what DCT is? – numbermaniac May 13 '17 at 0:32
• The Wikipedia page shows several different versions. Can you please clarify? Voting to close for now. – Nick Clifford May 13 '17 at 1:04
• @Nick Clifford lets assume DCT-II from that wikipedia page, as it clearly says: The DCT-II is probably the most commonly used form, and is often simply referred to as "the DCT" – user69099 May 13 '17 at 1:35
• Please add some test cases. – Shaggy May 13 '17 at 2:33
• @xakepp35 how do we know if our program is correct, if we don't have at least two test cases to test it against? Specifically, as long as you produce the same results as what? – Stephen May 18 '17 at 1:13

# Mathematica, 79 bytes

Tr/@(s=#;x=Length@s;Table[s[[n+1]]Cos[Pi/x(n+1/2)#],{n,0,x-1}]&/@Range[0,x-1])&


Since there are (currently) no test cases, as it turns out, there is a built-in for this: FourierDCT. However, the function "normalizes" the result by dividing it by a value before returning it. Thankfully, the documentation specifies that this division factor is just the square root of the input list's length (for DCT-II).

So, we can verify our results by multiplying the output of FourierDCT by the square root of the length of our input list. For anyone else who tries this problem, here are some test cases:

FourierDCT@{1,3,2,4}*Sqrt@4

{10, -2.38896, 0, -2.07193}

FourierDCT@{1,1,1,1,1}*Sqrt@5

{5, 0, 0, 0, 0}

Note that I've slightly edited the output of the last one; the last four values output as very small decimals, around 10-16, due to floating point precision errors around 0. My answer does give 0 for them, if you convert it to a decimal.

This solution gives exact answers when Mathematica can do so. If that's not okay and you want decimals, let me know.

# Jelly, 15 bytes

LḶµ+.×þ÷L-*Ḟ³æ.


Try it online!

This uses the identity that cos(pi * x) = real(e^(i * pi * x)) = real((-1)^x)

## Explanation

LḶµ+.×þ÷L-*Ḟ³æ.  Input: array A
L                Length
Ḷ               Lowered range - [0, 1, ..., len(A)-1]
µ              Begin a new monadic chain
×þ          Outer product using multiplication
÷         Divide by
L        Length
-*      Compute (-1)^x
Ḟ     Real part
æ.  Dot product with
³    The first argument (A)


## Clojure, 100 bytes

#(for[c[(count %)]k(range c)](apply +(map(fn[n x](*(Math/cos(*(/ c)Math/PI(+ n 0.5)k))x))(range)%)))


Sometimes it is difficult to be creative.

# C (gcc), 88 83 78 77 bytes

k;d(a,b,c)float*a,*b;{for(k=c*c;k--;)b[k/c]+=a[k%c]*cpow(-1,k/c*(.5+k%c)/c);}


Try it online!

This assumes the output vector is cleared before use.

Thanks to @Peter Taylor for -5

# Octave, 49 46 45 42 bytes

Here we consider the DCT as a matrix multiplication. The matrix is basically the entry-wise cosine of a rank 1 matrix which is very simple to construct.

Thanks @Cowsquack for -1 byte, thanks @ceilingcat for another -3 bytes!

@(x)cos(pi/(n=numel(x))*(0:n-1)'*(.5:n))*x


Try it online!

• Can you use .5 instead of 1/2? – user41805 Aug 9 '17 at 11:21
• D'oh, of course! – flawr Aug 9 '17 at 11:25
• Oh of course, thanks a lot! – flawr Sep 20 '17 at 16:15
• Suggest rows() instead of numel() – ceilingcat Oct 19 '19 at 18:55

# Mathematica, 51 49 bytes

Re@Total[I^Array[(2##+#2)/n&,{n=Length@#,n},0]#]&


Try it online! (Using Mathics)

Based on my solution to DFT.