(* DCT *)
(* $Id: dct.ml,v 1.1 2003/01/12 15:33:00 berke Exp $ *)
(* By Berke Durak *)
open Unix
open Bigarray
let pi = 3.14159265358979323846
let coef_8x8 =
let c = Array2.create float64 c_layout 8 8 in
for i = 0 to 7 do
let s = if i = 0 then sqrt 0.125 else 0.5
in
for j = 0 to 7 do
c.{i,j} <- s *. cos ((pi /. 8.0) *. (float_of_int i) *. ((float_of_int j) +. 0.5))
done
done;
c
let clip x =
let y = x + 128 in
if y < 0 then 0
else if y > 255 then 255
else y
let dct_8x8 tmp block =
for i = 0 to 7 do
for j = 0 to 7 do
let s = ref 0.0 in
for k = 0 to 7 do
s := !s +. coef_8x8.{j,k} *. (float_of_int block.{i,k})
done;
tmp.{i,j} <- !s;
done
done;
for j = 0 to 7 do
for i = 0 to 7 do
let s = ref 0.0 in
for k = 0 to 7 do
s := !s +. coef_8x8.{i,k} *. tmp.{k,j}
done;
block.{i,j} <- clip (int_of_float (floor (!s +. 0.499999)));
done
done
let idct_8x8 tmp block =
for i = 0 to 7 do
for j = 0 to 7 do
let s = ref 0.0 in
for k = 0 to 7 do
s := !s +. coef_8x8.{j,k} *. (float_of_int (block.{i,k} - 128))
done;
tmp.{i,j} <- !s;
done
done;
for j = 0 to 7 do
for i = 0 to 7 do
let s = ref 0.0 in
for k = 0 to 7 do
s := !s +. coef_8x8.{i,k} *. tmp.{k,j}
done;
block.{i,j} <- clip (int_of_float (floor (!s +. 0.499999)));
done
done
let write_pgm a fn =
let m = Array2.dim1 a
and n = Array2.dim2 a in
let oc = open_out_bin fn in
Printf.fprintf oc "P5 %d %d 255\n" n m;
for i = 0 to m - 1 do
for j = 0 to n - 1 do
output_char oc (Char.chr a.{i,j})
done
done;
close_out oc
let blit s d i0 j0 i1 j1 m n =
for i = 0 to m - 1 do
for j = 0 to n - 1 do
d.{i1 + i,j1 + j} <- s.{i0 + i, j0 + j}
done
done
let dimensions a = (Array2.dim1 a, Array2.dim2 a)
let create = Array2.create int8_unsigned c_layout
let createf = Array2.create float64 c_layout
let shuffle a =
let (m,n) = dimensions a in
let b = create m n in
let p = m / 8
and q = n / 8
in
for i = 0 to m - 1 do
for j = 0 to n - 1 do
b.{(i mod 8) * p + i / 8,(j mod 8) * q + j / 8} <- a.{i,j}
done
done;
b
let proc m n fn =
let fd = openfile fn [O_RDONLY] 0 in
let a = Array2.map_file fd int8_unsigned c_layout false m n in
let b = create m n in
let p = m / 8
and q = n / 8
in
let temp = createf 8 8
and block = create 8 8
and c = create m n
and d = create m n
in
for i = 0 to p - 1 do
for j = 0 to q - 1 do
blit a block (i * 8) (j * 8) 0 0 8 8;
dct_8x8 temp block;
blit block c 0 0 (i * 8) (j * 8) 8 8;
idct_8x8 temp block;
blit block d 0 0 (i * 8) (j * 8) 8 8;
done
done;
write_pgm c "/tmp/alpha.pgm";
write_pgm d "/tmp/beta.pgm";
let d = shuffle c in
write_pgm d "/tmp/gamma.pgm"
let _ =
proc (int_of_string Sys.argv.(1)) (int_of_string Sys.argv.(2)) Sys.argv.(3)