{-# OPTIONS_HADDOCK hide #-}

{-# LANGUAGE NoImplicitPrelude #-}

module Imj.Geo.Discrete.Resample
( resampleWithExtremities
) where

import           Imj.Prelude

import           Data.List( length )

import           Imj.Util( replicateElements )

{- | Resamples a list, using the analogy where a list
is seen as a uniform sampling of a geometrical segment.

With a uniform sampling strategy, for an input of length $$n$$, and a desired
output of length $$m$$:

* /Regular/ samples are repeated $$r = \lfloor {m \over n} \rfloor$$ times.
* /Over-represented/ samples are repeated $$r + 1$$ times.

If $$m'$$ is the number of over-represented samples,

\begin{alignedat}{2} m &= r*n + m' \\ \implies \quad m' &= m - r*n \end{alignedat}

We can chose over-represented samples in at least two different ways:

* __Even spread__ :

* Given a partition of the input continuous interval $$[\,0, length]\,$$
in $$m'$$ equal-length intervals, the over-represented samples are located at
the (floored) centers of these intervals.

* More precisely, over-represented samples indexes are:

$\biggl\{ a + \Bigl\lfloor {1 \over 2} + { n-1-a \over m-1 } * s \Bigl\rfloor \mid s \in [\,0\,..\,m'-1] \;,\; a = {1 \over 2} * {n \over m'} \biggl\}$

* Example : for a length 5 input, and 2 over-represented samples:

@
input samples:   -----

over-represented samples:    - -
@

* __"Even with extremities" spread__:

* The first and last over-represented samples match
with an input extremity. The rest of the over-represented samples are positionned
"regularly" in-between the first and last. An exception is made when there is only one
over-represented sample : in that case it is placed in the middle.

* More precisely, over-represented samples indexes are:

$if \; m' == 1 : \biggl\{ \Bigl\lfloor {n-1 \over 2} \Bigl\rfloor \biggl\}$

$otherwise : \biggl\{ \Bigl\lfloor {1 \over 2} + {n-1 \over m'-1}*s \Bigl\rfloor \mid s \in [\,0,m'-1]\, \biggl\}$

* Example : for a length 5 input, and 2 over-represented samples:

@
input samples:   -----

over-represented samples:   -   -
@

/As its name suggests, this function uses the "even with extremities" spread./

/For clarity, the variable names used in the code match the ones in the documentation./
-}
resampleWithExtremities :: [a]
-- ^ Input
-> Int
-- ^ $$n$$ : input length. It is expected that $$0 <= n <=$$ @length input@
-> Int
-- ^ $$m$$ : output length. It is expected that $$0 <= m$$.
-> [a]
-- ^ Output :
--
-- * when $$m < n$$, it is a /downsampled/ version of the input,
-- * when $$m > n$$, it is an /upsampled/ version of the input.
resampleWithExtremities input n m
| assert (m >= 0) m == n = input
| otherwise =
let r = quot m n
m' = m - (r * n)
res
| m' == 0   = replicateElements r input
| otherwise = let overRepIdx = getOverRepIdx (assert (m' > 0) m') n 0
in  resampleRec m' n 0 (overRepIdx, 0) input r
in  assert (verifyResample input m res) res

resampleRec :: Int
-- ^ over-represented samples count
-> Int
-- ^ $$n$$ : input length.
-> Int
-- ^ current index
-> (Int, Int)
-- ^ (next overrepresentation index, count of over-represented samples sofar)
-> [a]
-- ^ the list to be resampled
-> Int
-- ^  $$r = floor(m/n)$$ : every sample will be replicated
--   $$r$$ times, or $$r + 1$$ times if distance to next overrepresentation == 0
-> [a]
resampleRec _ _ _ _ [] _ = []
resampleRec m' n curIdx (overRepIdx, s) l@(_:_) r =
let (nCopies, nextState)
-- This commented guard was used to debug cases where the assert on the line after would fail
--        | overIdx < curIdx = error ("\noverIdx " ++ show overIdx ++ "\ncurIdx  " ++ show curIdx ++ "\nm' " ++ show m' ++ "\nn " ++ show n ++ "\ns " ++ show s)
| assert (overRepIdx >= curIdx) overRepIdx == curIdx
= let nextS = succ s
nextOverRepIdx = getOverRepIdx m' n nextS
in  (succ r, (nextOverRepIdx, nextS))
| otherwise = (r     , (overRepIdx    , s))
in replicate nCopies (head l) ++ resampleRec m' n (succ curIdx) nextState (tail l) r

-- | Returns maxBound when there is no over-representation
getOverRepIdx :: Int -> Int -> Int -> Int
getOverRepIdx m' n s
| m' >  1 = floor( 0.5 + (fromIntegral ((n - 1) * s) :: Float) / fromIntegral (m'-1))
| m' == 1 = if s == 0
then
quot n 2
else
maxBound
| otherwise = assert (m' == 0) maxBound

verifyResample :: [a]
-- ^ the input
-> Int
-- ^ the number of samples
-> [a]
-- ^ the output
-> Bool
verifyResample input nSamples resampled
| nSamples == length resampled = True
| otherwise                    = error \$ "\ninput    " ++ show (length input) ++
"\nnSamples " ++ show nSamples ++
"\nactual   " ++ show (length resampled)