module type T_Embed = sig
  type embed = float * float
  val string_of_embed : embed -> string
  val emb0 : embed
  val embed : float -> float -> embed 
  val emb_x : embed -> float
  val emb_y : embed -> float
  val embed_int : int -> int -> embed 
  val emb_x_int : embed -> int
  val emb_y_int : embed -> int
  val square : float -> float
  val scalaire : embed -> embed -> float
  val norm2 : embed -> embed -> float
  val norm : embed -> embed -> float
  val is_in_dist : embed -> embed -> float -> bool
  val middle : embed -> embed -> embed
  val move_alpha_1 : embed -> embed -> embed
  val move_alpha_2 : embed -> embed -> embed
  val alpha_0 : float -> embed -> embed -> embed
  val is_on : float -> float -> embed -> embed -> embed option -> float option
  val pi : float
  val degree_of_rad : float -> float
  val rad_of_degree : float -> float
  module Matrix : sig
    type vector = float array
    type matrix = float array array
    val vector : float -> float -> float -> vector
    val matrix : vector -> vector -> vector -> matrix
    val string_of_matrix : matrix -> string
    val vector_value : vector -> int -> float
    val matrix_line : matrix -> int -> vector
    val matrix_column : matrix -> int -> vector
    val matrix_transpose : matrix -> matrix
    val mult_add_2_vectors : vector -> vector -> float
    val mult_matrix_by_vector : matrix -> vector -> vector
    val mult_2_matrix : matrix -> matrix -> matrix
    val sample_vector : vector
    val sample_matrix : matrix
    val vector_from_embed : embed -> vector
    val vector_to_embed : vector -> embed
    val identity : matrix
    val translation : embed -> matrix
    val rotation_origin : float -> matrix
    val rotation : embed -> float -> matrix
    val homothetie_origin : float -> float -> matrix
    val homothetie : embed -> float -> float -> matrix
    val compose : matrix -> matrix -> matrix
    val apply : matrix -> embed -> embed
  end
end
module Embed : T_Embed = struct
  type embed = float * float
  let string_of_embed (x,y) = Printf.sprintf "<%f,%f>" x y
  let emb0 = (0.0,0.0)
  let embed x y = (x,y)
  let emb_x (x,y) = x
  let emb_y (x,y) = y
  let embed_int x y = (float_of_int x,float_of_int y)
  let emb_x_int (x,y) = int_of_float x
  let emb_y_int (x,y) = int_of_float y
  let square x = x *. x
  let scalaire (x,y) (x',y') = x *. x' +. y *. y'
  let scalaireT (x,y) (x',y') = x *. y' -. x' *. y
  let norm2 (x,y) (x',y') = square (x-.x') +. square (y-.y')
  let norm e e' = sqrt (norm2 e e')
  let is_in_dist e e' n = norm2 e e' < square n
  let middle (x,y) (x',y') = ((x+.x')/.2., (y+.y')/.2.)
  let move_alpha_1 (_,_) (x,y) = x,y
  let move_alpha_2 (_,_) (x,y) = x,y
  let alpha_0 d (x,y) (mx,my) =
    let dab = norm (x,y) (mx,my) in
    (mx+.(y-.my)*.d/.dab, my+.(mx-.x)*.d/.dab)
  let minmax seldist x y =
    if x>y then (y-.seldist, x+.seldist)
    else (x-.seldist, y+.seldist)
  let is_on min_dist sel_dist (x,y) (xa,ya) = function 
    | None ->
      let dist = norm2 (x,y) (xa,ya) in
      if dist < square min_dist then Some dist else None
    | Some (xb,yb) -> 
      let (xmin, xmax) = minmax sel_dist xa xb
      and (ymin, ymax) = minmax sel_dist ya yb in
      if x<xmin or x>xmax or y<ymin or y>ymax
      then None
      else let dx = (xb-.xa) and dy = (yb-.ya) in
      let dab = norm2 (xa,ya) (xb,yb)
      and f x y = x*.x/.y in
      let distT = f (scalaireT (x-.xa, y-.ya) (dx, dy)) dab
      and dist = f (scalaire (x-.xa, y-.ya) (dx, dy)) dab
      in if distT < square min_dist && dist < dab then Some dist else None
  let pi = 2. *. (acos 0.)
  let degree_of_rad a =
    let d = a *. 180. /. pi
    in if d < 0. then 360. +. d else d
  let rad_of_degree a = a *. pi /. 180.

  module Matrix = struct
    type vector = float array
    type matrix = float array array
    let vector x y z = [|x; y; z|]
    let matrix x y z = [|x; y; z|]
    let string_of_matrix m = 
      Printf.sprintf "{%f, %f, %f; %f, %f, %f; %f, %f, %f}"
        m.(0).(0) m.(0).(1) m.(0).(2)
        m.(1).(0) m.(1).(1) m.(1).(2)
        m.(2).(0) m.(2).(1) m.(2).(2)
    let vector_value v n = v.(n-1)
    let matrix_line m n = m.(n-1)
    let matrix_column m n =
      vector m.(0).(n-1) m.(1).(n-1) m.(2).(n-1)
    let matrix_transpose m = matrix
      (matrix_column m 1) (matrix_column m 2) (matrix_column m 3)
    let mult_add_2_vectors x y =
      x.(0) *. y.(0) +. x.(1) *. y.(1) +. x.(2) *. y.(2)
    let mult_matrix_by_vector m v = matrix
      (mult_add_2_vectors (matrix_line m 1) v)
      (mult_add_2_vectors (matrix_line m 2) v)
      (mult_add_2_vectors (matrix_line m 3) v)
    let mult_2_matrix m n = matrix_transpose (matrix
      (mult_matrix_by_vector m (matrix_column n 1))
      (mult_matrix_by_vector m (matrix_column n 2))
      (mult_matrix_by_vector m (matrix_column n 3)))
    let sample_vector = vector 1. 2. 3.
    let sample_matrix = matrix
      (vector 1. 2. 3.) (vector 4. 5. 6.) (vector 7. 8. 9.)
    let vector_from_embed (x, y) = vector x y 1.0
    let vector_to_embed v =
      (vector_value v 1) /. (vector_value v 3),
      (vector_value v 2) /. (vector_value v 3)
    let identity = matrix
      (vector 1.0 0.0 0.0)
      (vector 0.0 1.0 0.0)
      (vector 0.0 0.0 1.0)
    let translation (a,b) =
      Printf.printf "translation of %f, %f.\n" a b; flush stdout; matrix
      (vector 1.0 0.0 a)
      (vector 0.0 1.0 b)
      (vector 0.0 0.0 1.0)
    let rotation_origin a =
      let a = rad_of_degree a in matrix
      (vector (cos a) (-1.0 *. (sin a)) 0.0)
      (vector (sin a) (cos a) 0.0)
      (vector 0.0 0.0 1.0)
    let rotation (x,y) a = mult_2_matrix (translation (-.x, -.y))
        (mult_2_matrix (rotation_origin a) (translation (x, y)))
    let homothetie_origin a b = matrix
      (vector a 0.0 0.0)
      (vector 0.0 b 0.0)
      (vector 0.0 0.0 1.0)
    let homothetie (x,y) a b = mult_2_matrix (translation (-.x, -.y))
        (mult_2_matrix (homothetie_origin a b) (translation (x, y)))
    let compose m n = mult_2_matrix m n
    let apply m (x,y) = vector_to_embed
      (mult_matrix_by_vector m (vector_from_embed (x, y)))
  end
end