Simple OCaml implementation of the Discrete Fourier Transform for complex-valued signals.

dft.ml

(* complex num ops*)
type complex = {re: float; im: float}

let complex_add a b = {re = a.re +. b.re; im = a.im +. b.im} 
let complex_sub a b = {re = a.re -. b.re; im = a.im -. b.im}
let complex_mul a b = {re = a.re *. b.re -. a.im *. b.im; im = a.re *. b.im +. a.im *. b.re}
(* exponential of a complex number given angle in radians*)
let complex_exp theta = {re = cos theta; im = sin theta}

(* discrete fourier transform*)
let dft x = 
  let n = Array.length x in
  let x_k k =
    let sum = ref {re = 0.0; im = 0.0} in
    for n_i = 0 to n - 1 do
      let e = complex_exp (-. 2.0 *. Float.pi *. float_of_int k *. float_of_int n_i /. float_of_int n) in
      sum := complex_add !sum (complex_mul x.(n_i) e)
    done;
    !sum 
  in 
  Array.init n x_k 

(* example use *)
let () =
  let signal = [|{re = 1.0; im = 0.0}; {re = 2.0; im = 0.0}; {re = 3.0; im = 0.0}; {re = 4.0; im = 0.0}|] in
  let result = dft signal in 
  Array.iter (fun c -> Printf.printf "(%f, %f)\n" c.re c.im) result