-- | <http://paulbourke.net/>
module Sound.SC3.Data.Math.Bourke where

-- | 2-element /h/ transform.
h_transform_2 :: Num t => t -> ((t,t) -> (t,t)) -> (t,t) -> (t,t)
h_transform_2 h f (x,y) =
    let (x',y') = f (x,y)
    in (x + h * x',y + h * y')

-- | 3-element /h/ transform.
h_transform_3 :: Num t => t -> ((t,t,t) -> (t,t,t)) -> (t,t,t) -> (t,t,t)
h_transform_3 h f (x,y,z) =
    let (x',y',z') = f (x,y,z)
    in (x + h * x',y + h * y',z + h * z')

-- | <http://paulbourke.net/fractals/lorenz/>
lorenz :: Num t => t -> t -> t -> (t,t,t) -> (t,t,t)
lorenz a b c (x,y,z) =
    (a * (y - x)
    ,x * (b - z) - y
    ,x * y - c * z)

{- | <http://paulbourke.net/fractals/lorenz/>

> import Sound.SC3.Plot {- hsc3-plot -}

> let l = iterate (lorenz_h 0.01 10 28 (8/3)) (0.1,0.0,0.0)
> in plot_p3_ln [take 5000 l]

> let {l = iterate (lorenz_h 0.01 10 28 (8/3)) (0.1,0.0,0.0)
>     ;f (x,_,z) = (x,z)}
> in plot_p2_ln [take 15000 (map f l)]

-}
lorenz_h :: Num t => t -> t -> t -> t -> (t,t,t) -> (t,t,t)
lorenz_h h a b c = h_transform_3 h (lorenz a b c)

-- | <http://paulbourke.net/fractals/rossler/>
rossler :: Num t => t -> t -> t -> (t, t, t) -> (t, t, t)
rossler a b c (x,y,z) =
    (negate y - z
    ,x + a * y
    ,b + z * (x - c))

{- | <http://paulbourke.net/fractals/rossler/>

> plot_p3_ln [take 5000 (iterate (rossler_h 0.02 0.2 0.2 5.7) (0.1,0,0))]

-}
rossler_h :: Num t => t -> t -> t -> t -> (t,t,t) -> (t,t,t)
rossler_h h a b c = h_transform_3 h (rossler a b c)

{- | <http://paulbourke.net/fractals/peterdejong/>

> let pdj a b c d =
>   let vw x = plot_p2_pt [take 15000 x]
>   in vw (iterate (peter_de_jong a b c d) (-0.72,-0.64))

> pdj 1.4 (-2.3) 2.4 (-2.1)
> pdj 2.01 (-2.53) 1.61 (-0.33)
> pdj (-2.7) (-0.09) (-0.86) (-2.2)
> pdj (-2.24) 0.43 (-0.65) (-2.43)
> pdj (-2.0) (-2.0) (-1.2) 2.0

-}
peter_de_jong :: Floating t => t -> t -> t -> t -> (t, t) -> (t, t)
peter_de_jong a b c d (x,y) =
    (sin (a * y) - cos (b * x)
    ,sin (c * x) - cos (d * y))

{- | <http://paulbourke.net/fractals/clifford/>

> let clf a b c d =
>   let vw x = plot_p2_pt [take 12500 x]
>   in vw (iterate (clifford a b c d) (-0.72,-0.64))

> clf (-1.4) (1.6) (1.0) (0.7)
> clf (1.1) (-1.0) (1.0) (1.5) {- not as pb indicates -}
> clf (1.6) (-0.6) (-1.2) (1.6)
> clf (1.7) (1.7) (0.6) (1.2)

-}
clifford :: Floating t => t -> t -> t -> t -> (t, t) -> (t, t)
clifford a b c d (x,y) =
    (sin (a * y) + c * cos (a * x)
    ,sin (b * x) + d * cos (b * y))