module Geometry.Trajectory.PathSimplifier.RamerDouglasPeucker (
      simplifyTrajectoryRdp
    , simplifyTrajectoryRdpBy
) where



import           Data.Ord
import           Data.Sequential
import           Data.Vector     (Vector)
import qualified Data.Vector     as V

import Geometry.Core



-- | Simplify a path by dropping unnecessary points, using the the
-- [Ramer-Douglas-Peucker algorithm]
-- (https://en.wikipedia.org/wiki/Ramer%E2%80%93Douglas%E2%80%93Peucker_algorithm).
-- The larger the tolerance distance, the simpler the result will be.
--
-- This is very useful in conjunction with 'bezierSmoothen': first drop the
-- redundancies, then smoothen using Bezier curves again, to yield a result
-- visually similar to the original data, but with a much smaller data footprint
-- (SVGs can become huge!).
--
-- If your trajectory contains more than just the points you want to simplify on,
-- use 'simplifyTrajectoryRdpBy'.
--
-- <<docs/interpolation/3_simplify_path_rdp.svg>>
simplifyTrajectoryRdp
    :: Sequential vector
    => Double      -- ^ Discard points closer than \(\varepsilon\) to the connecting line of two points. Larger values yield simpler results.
    -> vector Vec2 -- ^ Trajectory
    -> Vector Vec2 -- ^ Simplified trajectory
simplifyTrajectoryRdp :: forall (vector :: * -> *).
Sequential vector =>
Double -> vector Vec2 -> Vector Vec2
simplifyTrajectoryRdp Double
epsilon = Double -> (Vec2 -> Vec2) -> vector Vec2 -> Vector Vec2
forall (vector :: * -> *) a.
Sequential vector =>
Double -> (a -> Vec2) -> vector a -> Vector a
simplifyTrajectoryRdpBy Double
epsilon Vec2 -> Vec2
forall a. a -> a
id

-- | 'simplifyTrajectoryRdp', but allows specifying a function for how to extract the points
-- to base simplifying on from the input.
--
-- This is useful when your trajectory contains metadata, such as the velocity at
-- each point.
simplifyTrajectoryRdpBy
    :: Sequential vector
    => Double      -- ^ Discard points closer than \(\varepsilon\) to the connecting line of two points. Larger values yield simpler results.
    -> (a -> Vec2) -- ^ Extract the relevant 'Vec2' to simplify on
    -> vector a    -- ^ Trajectory
    -> Vector a    -- ^ Simplified trajectory
simplifyTrajectoryRdpBy :: forall (vector :: * -> *) a.
Sequential vector =>
Double -> (a -> Vec2) -> vector a -> Vector a
simplifyTrajectoryRdpBy Double
epsilon a -> Vec2
vec2in = Vector a -> Vector a
go (Vector a -> Vector a)
-> (vector a -> Vector a) -> vector a -> Vector a
forall b c a. (b -> c) -> (a -> b) -> a -> c
. vector a -> Vector a
forall a. vector a -> Vector a
forall (f :: * -> *) a. Sequential f => f a -> Vector a
toVector
  where
    go :: Vector a -> Vector a
go Vector a
points | Vector a -> Int
forall a. Vector a -> Int
V.length Vector a
points Int -> Int -> Bool
forall a. Ord a => a -> a -> Bool
<= Int
2 = Vector a
points
    go Vector a
points
      = let start :: a
start = Vector a -> a
forall a. Vector a -> a
V.head Vector a
points
            end :: a
end = Vector a -> a
forall a. Vector a -> a
V.last Vector a
points
            isCyclic :: Bool
isCyclic = a -> Vec2
vec2in a
start Vec2 -> Vec2 -> Bool
forall a. Eq a => a -> a -> Bool
== a -> Vec2
vec2in a
end

            (Double
dMax, Int
iMax)
                | Bool
isCyclic = ((Double, Int) -> (Double, Int) -> Ordering)
-> Vector (Double, Int) -> (Double, Int)
forall a. (a -> a -> Ordering) -> Vector a -> a
V.maximumBy (((Double, Int) -> Double)
-> (Double, Int) -> (Double, Int) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Double, Int) -> Double
forall a b. (a, b) -> a
fst) ((Int -> a -> (Double, Int)) -> Vector a -> Vector (Double, Int)
forall a b. (Int -> a -> b) -> Vector a -> Vector b
V.imap (\Int
i a
p -> (Vec2 -> Double
norm (a -> Vec2
vec2in a
p Vec2 -> Vec2 -> Vec2
forall v. VectorSpace v => v -> v -> v
-. a -> Vec2
vec2in a
start), Int
i)) Vector a
points)
                | Bool
otherwise =
                    let connectingLine :: Line
connectingLine = Vec2 -> Vec2 -> Line
Line (a -> Vec2
vec2in a
start) (a -> Vec2
vec2in a
end)
                    in ((Double, Int) -> (Double, Int) -> Ordering)
-> Vector (Double, Int) -> (Double, Int)
forall a. (a -> a -> Ordering) -> Vector a -> a
V.maximumBy (((Double, Int) -> Double)
-> (Double, Int) -> (Double, Int) -> Ordering
forall a b. Ord a => (b -> a) -> b -> b -> Ordering
comparing (Double, Int) -> Double
forall a b. (a, b) -> a
fst) ((Int -> a -> (Double, Int)) -> Vector a -> Vector (Double, Int)
forall a b. (Int -> a -> b) -> Vector a -> Vector b
V.imap (\Int
i a
p -> (Vec2 -> Line -> Double
distanceFromLine (a -> Vec2
vec2in a
p) Line
connectingLine, Int
i)) Vector a
points)
        in if Double
dMax Double -> Double -> Bool
forall a. Ord a => a -> a -> Bool
<= Double
epsilon
            -- Points are all too close: remove them
            then [a] -> Vector a
forall a. [a] -> Vector a
V.fromList [a
start, a
end]
            -- Points far enough away: recurse
            else let before :: Vector a
before = Int -> Vector a -> Vector a
forall a. Int -> Vector a -> Vector a
V.take (Int
iMaxInt -> Int -> Int
forall a. Num a => a -> a -> a
+Int
1) Vector a
points -- +1 so this includes the iMax point
                     after :: Vector a
after = Int -> Vector a -> Vector a
forall a. Int -> Vector a -> Vector a
V.drop Int
iMax Vector a
points
                 in Vector a -> Vector a
go Vector a
before Vector a -> Vector a -> Vector a
forall a. Semigroup a => a -> a -> a
<> Int -> Vector a -> Vector a
forall a. Int -> Vector a -> Vector a
V.drop Int
1 (Vector a -> Vector a
go Vector a
after) -- Don’t add the middle twice