{-# LINE 1 "Physics/Hipmunk/Shape.hsc" #-}
-----------------------------------------------------------------------------
{-# LINE 2 "Physics/Hipmunk/Shape.hsc" #-}
-- |
-- Module      :  Physics/Hipmunk/Shape.hsc
-- Copyright   :  (c) Felipe A. Lessa 2008
-- License     :  MIT (see LICENSE)
--
-- Maintainer  :  felipe.lessa@gmail.com
-- Stability   :  provisional
-- Portability :  portable (needs FFI)
--
-- Shapes used for collisions, their properties and some useful
-- polygon functions.
--
-----------------------------------------------------------------------------

module Physics.Hipmunk.Shape
    (-- * Shapes
     Shape,
     ShapeType(..),
     newShape,

     -- * Properties
     -- ** Collision type
     CollisionType,
     getCollisionType,
     setCollisionType,
     -- ** Group
     Group,
     getGroup,
     setGroup,
     -- ** Layers
     Layers,
     getLayers,
     setLayers,
     -- ** Elasticity
     Elasticity,
     getElasticity,
     setElasticity,
     -- ** Friction
     Friction,
     getFriction,
     setFriction,
     -- ** Surface velocity
     SurfaceVel,
     getSurfaceVel,
     setSurfaceVel,

     -- * Utilities
     getBody,
     momentForCircle,
     momentForPoly,
     shapeQuery,

     -- ** For polygons
     -- $polygon_util
     Segment,
     Intersection(..),
     epsilon,
     (.==.),
     isLeft,
     isClockwise,
     isConvex,
     intersects,
     polyReduce,
     polyCenter,
     convexHull
    )
    where

import Data.List (foldl', sortBy)
import Foreign hiding (rotate, new)
import Foreign.C

{-# LINE 74 "Physics/Hipmunk/Shape.hsc" #-}

import Physics.Hipmunk.Common
import Physics.Hipmunk.Internal

-- | There are three types of shapes that can be attached
--   to bodies:
data ShapeType =
    -- | A circle is the fastest collision type. It also
    --   rolls smoothly.
    Circle {radius :: !CpFloat}

    -- | A line segment is meant to be used as a static
    --   shape. (It can be used with moving bodies, however
    --   two line segments never generate collisions between
    --   each other.)
  | LineSegment {start     :: !Position,
                 end       :: !Position,
                 thickness :: !CpFloat}

    -- | Polygons are the slowest of all shapes but
    --   the most flexible. The list of vertices must form
    --   a convex hull with clockwise winding.
    --   Note that if you want a non-convex polygon you may
    --   add several convex polygons to the body.
  | Polygon {vertices :: ![Position]}
    deriving (Eq, Ord, Show)


-- | @newShape b type off@ creates a new shape attached to
--   body @b@ at offset @off@. Note that you have to
--   add the shape to a space otherwise it won't generate
--   collisions.
newShape :: Body -> ShapeType -> Position -> IO Shape
newShape body@(B b) (Circle r) off =
  withForeignPtr b $ \b_ptr ->
  with off $ \off_ptr ->
  mallocForeignPtrBytes (80) >>= \shape ->
{-# LINE 111 "Physics/Hipmunk/Shape.hsc" #-}
  withForeignPtr shape $ \shape_ptr -> do
    wrCircleShapeInit shape_ptr b_ptr off_ptr r
    return (S shape body)

newShape body@(B b) (LineSegment p1 p2 r) off =
  withForeignPtr b $ \b_ptr ->
  with (p1+off) $ \p1off_ptr ->
  with (p2+off) $ \p2off_ptr ->
  mallocForeignPtrBytes (112) >>= \shape ->
{-# LINE 120 "Physics/Hipmunk/Shape.hsc" #-}
  withForeignPtr shape $ \shape_ptr -> do
    wrSegmentShapeInit shape_ptr b_ptr p1off_ptr p2off_ptr r
    return (S shape body)

newShape body@(B b) (Polygon verts) off =
  withForeignPtr b $ \b_ptr ->
  with off $ \off_ptr ->
  withArrayLen verts $ \verts_len verts_ptr ->
  mallocForeignPtrBytes (80) >>= \shape ->
{-# LINE 129 "Physics/Hipmunk/Shape.hsc" #-}
  withForeignPtr shape $ \shape_ptr -> do
    let verts_len' = fromIntegral verts_len
    wrPolyShapeInit shape_ptr b_ptr verts_len' verts_ptr off_ptr
    addForeignPtrFinalizer cpShapeDestroy shape
    return (S shape body)

foreign import ccall unsafe "wrapper.h"
    wrCircleShapeInit :: ShapePtr -> BodyPtr -> VectorPtr
                      -> CpFloat -> IO ()
foreign import ccall unsafe "wrapper.h"
    wrSegmentShapeInit :: ShapePtr -> BodyPtr -> VectorPtr
                       -> VectorPtr -> CpFloat -> IO ()
foreign import ccall unsafe "wrapper.h"
    wrPolyShapeInit :: ShapePtr -> BodyPtr -> CInt -> VectorPtr
                    -> VectorPtr -> IO ()
foreign import ccall unsafe "wrapper.h &cpShapeDestroy"
    cpShapeDestroy :: FunPtr (ShapePtr -> IO ())


-- | @getBody s@ is the body that this shape is associated
--   to. Useful especially in 'Physics.Hipmunk.Space.Callback'.
getBody :: Shape -> Body
getBody (S _ b) = b


-- | The collision type is used to determine which collision
--   'Physics.Hipmunk.Space.Callback' will be called. Its
--   actual value doesn't have a meaning for Chipmunk other
--   than the correspondence between shapes and the collision
--   pair functions you add. (default is zero)
type CollisionType = Word32
{-# LINE 160 "Physics/Hipmunk/Shape.hsc" #-}
getCollisionType :: Shape -> IO CollisionType
getCollisionType (S shape _) =
  withForeignPtr shape (\hsc_ptr -> peekByteOff hsc_ptr 24)
{-# LINE 163 "Physics/Hipmunk/Shape.hsc" #-}
setCollisionType :: Shape -> CollisionType -> IO ()
setCollisionType (S shape _) col =
  withForeignPtr shape $ \shape_ptr -> do
    (\hsc_ptr -> pokeByteOff hsc_ptr 24) shape_ptr col
{-# LINE 167 "Physics/Hipmunk/Shape.hsc" #-}

-- | Groups are used to filter collisions between shapes. If
--   the group is zero, then it imposes no restriction
--   to the collisions. However, if the group is non-zero then
--   the shape will not collide with other shapes in the same
--   non-zero group. (default is zero)
--
--   This is primarely used to create multi-body, multi-shape
--   objects such as ragdolls. It may be thought as a lightweight
--   alternative to creating a callback that filters the
--   collisions.
type Group = Word32
{-# LINE 179 "Physics/Hipmunk/Shape.hsc" #-}
getGroup :: Shape -> IO Group
getGroup (S shape _) =
  withForeignPtr shape (\hsc_ptr -> peekByteOff hsc_ptr 28)
{-# LINE 182 "Physics/Hipmunk/Shape.hsc" #-}
setGroup :: Shape -> Group -> IO ()
setGroup (S shape _) gr =
  withForeignPtr shape $ \shape_ptr -> do
    (\hsc_ptr -> pokeByteOff hsc_ptr 28) shape_ptr gr
{-# LINE 186 "Physics/Hipmunk/Shape.hsc" #-}

-- | Layers are similar to groups, but use a bitmask. For a collision
--   to occur, two shapes must have at least one layer in common.
--   In other words, @layer1 .&. layer2@ should be non-zero.
--   (default is @0xFFFF@)
--
--   Note that although this type may have more than 32 bits,
--   for portability you should only rely on the lower 32 bits.
type Layers = Word32
{-# LINE 195 "Physics/Hipmunk/Shape.hsc" #-}
getLayers :: Shape -> IO Layers
getLayers (S shape _) =
  withForeignPtr shape (\hsc_ptr -> peekByteOff hsc_ptr 32)
{-# LINE 198 "Physics/Hipmunk/Shape.hsc" #-}
setLayers :: Shape -> Layers -> IO ()
setLayers (S shape _) lay =
  withForeignPtr shape $ \shape_ptr -> do
    (\hsc_ptr -> pokeByteOff hsc_ptr 32) shape_ptr lay
{-# LINE 202 "Physics/Hipmunk/Shape.hsc" #-}

-- | The elasticity of the shape is such that @0.0@ gives no bounce
--   while @1.0@ give a \"perfect\" bounce. Note that due to
--   inaccuracies using @1.0@ or greater is not recommended.
--
--   The amount of elasticity applied during a collision is
--   calculated by multiplying the elasticity of both shapes.
--   (default is zero)
--
--   /IMPORTANT:/ by default no elastic iterations are done
--   when the space 'Physics.Hipmunk.Space.step's. This means
--   that all shapes react as they had zero elasticity.
--   So, if you want some elasticity, remember to call
--   'Physics.Hipmunk.Space.setElasticIterations' to something
--   greater than zero, maybe @10@.
type Elasticity = CpFloat
getElasticity :: Shape -> IO Elasticity
getElasticity (S shape _) =
  withForeignPtr shape (\hsc_ptr -> peekByteOff hsc_ptr 44)
{-# LINE 221 "Physics/Hipmunk/Shape.hsc" #-}
setElasticity :: Shape -> Elasticity -> IO ()
setElasticity (S shape _) e =
  withForeignPtr shape $ \shape_ptr -> do
    (\hsc_ptr -> pokeByteOff hsc_ptr 44) shape_ptr e
{-# LINE 225 "Physics/Hipmunk/Shape.hsc" #-}

-- | The friction coefficient of the shape according
--   to Coulumb friction model (i.e. @0.0@ is frictionless,
--   iron on iron is around @1.0@, and it could be greater
--   then @1.0@).
--
--   The amount of friction applied during a collision is
--   determined by multiplying the friction coefficient
--   of both shapes. (default is zero)
type Friction = CpFloat
getFriction :: Shape -> IO Friction
getFriction (S shape _) =
  withForeignPtr shape (\hsc_ptr -> peekByteOff hsc_ptr 48)
{-# LINE 238 "Physics/Hipmunk/Shape.hsc" #-}
setFriction :: Shape -> Friction -> IO ()
setFriction (S shape _) u =
  withForeignPtr shape $ \shape_ptr -> do
    (\hsc_ptr -> pokeByteOff hsc_ptr 48) shape_ptr u
{-# LINE 242 "Physics/Hipmunk/Shape.hsc" #-}

-- | The surface velocity of the shape. Useful to create
--   conveyor belts and players that move around. This
--   value is only used when calculating friction, not
--   collision. (default is zero)
type SurfaceVel = Vector
getSurfaceVel :: Shape -> IO SurfaceVel
getSurfaceVel (S shape _) =
  withForeignPtr shape (\hsc_ptr -> peekByteOff hsc_ptr 52)
{-# LINE 251 "Physics/Hipmunk/Shape.hsc" #-}
setSurfaceVel :: Shape -> SurfaceVel -> IO ()
setSurfaceVel (S shape _) sv =
  withForeignPtr shape $ \shape_ptr -> do
    (\hsc_ptr -> pokeByteOff hsc_ptr 52) shape_ptr sv
{-# LINE 255 "Physics/Hipmunk/Shape.hsc" #-}




-- | @momentForCircle m (ri,ro) off@ is the moment of inertia
--   of a circle of @m@ mass, inner radius of @ri@, outer radius
--   of @ro@ and at an offset @off@ from the center of the body.
momentForCircle :: CpFloat -> (CpFloat, CpFloat) -> Position -> CpFloat
momentForCircle m (ri,ro) off = (m/2)*(ri*ri + ro*ro) + m*(off `dot` off)
-- We recoded the C function to avoid FFI and unsafePerformIO
-- on this simple function.


-- | @momentForPoly m verts off@ is the moment of inertia of a
--   polygon of @m@ mass, at offset @off@ from the center of
--   the body and comprised of @verts@ vertices. This is similar
--   to 'shapePoly' (and the same restrictions for the vertices
--   apply as well).
momentForPoly :: CpFloat -> [Position] -> Position -> CpFloat
momentForPoly m verts off = (m*sum1)/(6*sum2)
  where
    verts' = if off /= 0 then map (+off) verts else verts
    (sum1,sum2) = calc (pairs (,) verts') 0 0

    calc a b c | a `seq` b `seq` c `seq` False = undefined
    calc []           acc1 acc2 = (acc1, acc2)
    calc ((v1,v2):vs) acc1 acc2 =
      let a = v2 `cross` v1
          b = v1 `dot` v1 + v1 `dot` v2 + v2 `dot` v2
      in calc vs (acc1 + a*b) (acc2 + a)
-- We recoded the C function to avoid FFI, unsafePerformIO
-- and a bunch of malloc + poke. Is it worth?

-- | Internal. For @l = [x1,x2,...,xn]@, @pairs f l@ is
--   @[f x1 x2, f x2 x3, ...,f xn x1]@.
pairs :: (a -> a -> b) -> [a] -> [b]
pairs f l = zipWith f l (tail $ cycle l)

-- | @shapeQuery shape p@ returns @True@ iff the point in
--   position @p@ (in world's coordinates) lies within
--   the shape @shape@.
shapeQuery :: Shape -> Position -> IO Bool
shapeQuery (S shape _) p =
  withForeignPtr shape $ \shape_ptr ->
  with p $ \p_ptr -> do
    i <- wrShapePointQuery shape_ptr p_ptr
    return (i /= 0)

foreign import ccall unsafe "wrapper.h"
    wrShapePointQuery :: ShapePtr -> VectorPtr -> IO CInt



-- $polygon_util
--   This section is inspired by @pymunk.util@,
--   a Python module made from <http://code.google.com/p/pymunk/>,
--   although implementations are quite different.
--
--   Also, unless noted otherwise all polygons are
--   assumed to be simple (i.e. no overlapping edges).

-- | The epsilon used in the algorithms below when necessary
--   to compare floats for \"equality\".
epsilon :: CpFloat
epsilon = 1e-25

-- | \"Equality\" under 'epsilon'. That is, @a .==. b@
--   if @abs (a - b) <= epsilon@.
(.==.) :: CpFloat -> CpFloat -> Bool
a .==. b = abs (a - b) <= epsilon

-- | A line segment.
type Segment = (Position, Position)

-- | /O(n)/. @isClockwise verts@ is @True@ iff @verts@ form
--   a clockwise polygon.
isClockwise :: [Position] -> Bool
isClockwise = (<= 0) . foldl' (+) 0 . pairs cross

-- | @isLeft (p1,p2) vert@ is
--
--    * @LT@ if @vert@ is at the left of the line defined by @(p1,p2)@.
--
--    * @EQ@ if @vert@ is at the line @(p1,p2)@.
--
--    * @GT@ otherwise.
isLeft :: (Position, Position) -> Position -> Ordering
isLeft (p1,p2) vert = compare 0 $ (p1 - vert) `cross` (p2 - vert)

-- | /O(n)/. @isConvex verts@ is @True@ iff @vers@ form a convex
--   polygon.
isConvex :: [Position] -> Bool
isConvex = foldl1 (==) . map (0 <) . filter (0 /=) . pairs cross . pairs (-)
-- From http://apocalisp.wordpress.com/category/programming/haskell/page/2/

-- | /O(1)/. @intersects seg1 seg2@ is the intersection between
--   the two segments @seg1@ and @seg2@. See 'Intersection'.
intersects :: Segment -> Segment -> Intersection
intersects (a0,a1) (b0,b1) =
    let u                = a1 - a0
        v@(Vector vx vy) = b1 - b0
        w@(Vector wx wy) = a0 - b0
        d = u `cross` v
        parallel = d .==. 0

        -- Parallel case
        collinear = all (.==. 0) [u `cross` w, v `cross` w]
        a_is_point = u `dot` u .==. 0
        b_is_point = v `dot` v .==. 0
        (Vector w2x w2y) = a1 - b0
        (a_in_b, a_in_b') = if vx .==. 0
                             then swap (wy/vy, w2y/vy)
                             else swap (wx/vx, w2x/vx)
            where swap t@(x,y) | x < y     = t
                               | otherwise = (y,x)

        -- Non-parallel case
        sI = v `cross` w / d
        tI = u `cross` w / d

        -- Auxiliary functions
        inSegment p (c0,c1)
            | vertical  = test (gy p) (gy c0, gy c1)
            | otherwise = test (gx p) (gx c0, gx c1)
            where
              vertical = gx c0 .==. gx c1
              (gx, gy) = (\(Vector x _) -> x, \(Vector _ y) -> y)
              test q (d0,d1) = any (inside q) [(d0,d1), (d1,d0)]
        inside n (l,r) = l <= n && n <= r

    in if parallel
       then case (collinear, a_is_point, b_is_point) of
             (False, _, _) ->
                 -- Parallel and non-collinear
                 IntNowhere

             (_, False, False) ->
                 -- Both are parallel, collinear segments
                 case (a_in_b > 1 || a_in_b' < 0,
                       max a_in_b 0, min a_in_b' 1) of
                   (True, _, _) -> IntNowhere
                   (_, i0, i1)
                       | i0 .==. i1 -> IntPoint p0
                       | otherwise  -> IntSegmt (p0,p1)
                       where p0 = b0 + v `scale` i0
                             p1 = b0 + v `scale` i1

             (_, True, True) ->
                 -- Both are points
                 if len (b0-a0) .==. 0
                 then IntPoint a0 else IntNowhere

             _ ->
                 -- One is a point, another is a segment
                 let (point,segment)
                         | a_is_point = (a0, (b0,b1))
                         | otherwise  = (b0, (a0,a1))
                 in if inSegment point segment
                    then IntPoint point else IntNowhere

       else if all (\x -> inside x (0,1)) [sI, tI]
            then IntPoint (a0 + u `scale` sI) else IntNowhere

-- | A possible intersection between two segments.
data Intersection = IntNowhere         -- ^ Don't intercept.
                  | IntPoint !Position -- ^ Intercept in a point.
                  | IntSegmt !Segment  -- ^ Share a segment.
                    deriving (Eq, Ord, Show)


-- | /O(n)/. @polyReduce delta verts@ removes from @verts@ all
--   points that have less than @delta@ distance
--   in relation to the one preceding it.
--
--   Note that a very small polygon may be completely \"eaten\"
--   if all its vertices are within a @delta@ radius from the
--   first.
polyReduce :: CpFloat -> [Position] -> [Position]
polyReduce delta = go
    where
      go (p1:p2:ps) | len (p2-p1) < delta = go (p1:ps)
                    | otherwise           = p1 : go (p2:ps)
      go other = other

-- | /O(n)/. @polyCenter verts@ is the position in the center
--   of the polygon formed by @verts@.
polyCenter :: [Position] -> Position
polyCenter verts = foldl' (+) 0 verts `scale` s
    where s = recip $ toEnum $ length verts


-- | /O(n log n)/. @convexHull verts@ is the convex hull of the
--   polygon defined by @verts@. The vertices of the convex
--   hulls are given in clockwise winding. The polygon
--   doesn't have to be simple.
--
--   Implemented using Graham scan, see
--   <http://cgm.cs.mcgill.ca/~beezer/cs507/3coins.html>.
convexHull :: [Position] -> [Position]
convexHull verts =
  let (p0,ps) = takeMinimum verts
      (_:p1:points) = p0 : sortBy (isLeft . (,) p0) ps

      -- points is going counterclockwise now.
      -- In go we use 'hull' with the last added
      -- vertex as the head, so our result is clockwise.

      -- Remove right turns
      go hull@(h1:h2:hs) (q1:qs) =
          case (isLeft (h2,h1) q1, hs) of
            (LT,_) -> go (q1:hull) qs    -- Left turn
            (_,[]) -> go (q1:hull) qs    -- Maintain at least 2 points
            _      -> go (h2:hs) (q1:qs) -- Right turn or straight
      go hull [] = hull
      go _ _     = error "Physics.Hipmunk.Shape.convexHull: never get here"

  in go [p1,p0] points


-- | Internal. Works like minimum but also returns the
--   list without it. The order of the list may be changed.
--   We have @fst (takeMinimum xs) == minimum xs@ and
--   @sort (uncurry (:) $ takeMinimum xs) == sort xs@
takeMinimum :: Ord a => [a] -> (a, [a])
takeMinimum [] = error "Physics.Hipmunk.Shape.takeMinimum: empty list"
takeMinimum (x:xs) = go x [] xs
    where
      go min_ acc (y:ys) | y < min_  = go y (min_:acc) ys
                         | otherwise = go min_ (y:acc) ys
      go min_ acc [] = (min_, acc)




{-

-- | /O((n+k)log n)/ [where /k/ is the number of intersections].
--   @intersections segs@ is the list of all intersections
--   found between the list @segs@ of line segments. Each
--   intersection is represented as two integers meant
--   to be interpreted as two indexes of @segs@ (zero being
--   the first line segment).
--
--   It is a implementation of the Bentley-Ottmann algorithm (see
--   <http://geometryalgorithms.com/Archive/algorithm_0108/algorithm_0108.htm>
--   for example), however intersection points are \"returned\"
--   as soon as they are found. That is, the WHNF of @intersections segs@
--   only needs to calculate the necessary to find the first
--   intersection, so if you want to know only if there is a
--   an intersection or not then you only need /O(n log n)/ time.
--
--   (Note that the @segs@ does not need to be a polygon at all.)
intersections :: [Segment] -> [InterIndexes]
intersections = bentleyOttmann . zip [0..]

type InterIndexes = (Intersection, SegmentIndex, SegmentIndex)


--- Basic data types
type SegmentIndex = Int
type SegmentArray = Array SegmentIndex Segment
type IndSeg = (SegmentIndex, Segment)
type Neighbors = (SegmentIndex, SegmentIndex)
data Event = EvStart !Position !SegmentIndex
           | EvEnd   !Position !SegmentIndex
           | EvInter !Position !SegmentIndex !SegmentIndex
             deriving (Eq)

instance Ord Event where
    e1 `compare` e2 = case evPos e1 `compare` evPos e2 of
                        EQ -> evIdent e1 `compare` evIdent e2
                        ot -> ot
        where
          evIdent (EvStart _ i)   = (-2, i)
          evIdent (EvEnd _ i)     = (-1, i)
          evIdent (EvInter _ i j) = (i, j)

evPos :: Event -> Position
evPos (EvStart p _)   = p
evPos (EvEnd p _)     = p
evPos (EvInter p _ _) = p

interErr :: a
interErr = error . ("Physics.Hipmunk.Shape.intersections: " ++)

interDebug :: Bool -> Bool
--interDebug = const False
interDebug = id

--- Event queue (a priority queue)
data EventQueue = EQNil
                | EQBranch !Event EventQueue EventQueue

eqSingle :: Event -> EventQueue
eqSingle ev = EQBranch ev EQNil EQNil

eqGet :: EventQueue -> (Event, EventQueue)
eqGet EQNil               = interErr "[eqGet] never get here"
eqGet (EQBranch ev q1 q2) = (ev, eqMerge q1 q2)

eqMerge :: EventQueue -> EventQueue -> EventQueue
eqMerge EQNil other = other
eqMerge other EQNil = other
eqMerge left@(EQBranch evL _ _) right@(EQBranch evR r1 r2)
    = case ev `compare` ev' of
        LT -> helper left right
        GT -> helper right left
        EQ -> eqMerge left (eqMerge r1 r2) -- Discard ev'!
    | ev <= ev' = helper left right
    | otherwise = helper right left
    where
      helper (EQBranch ev EQNil q2) r = EQBranch ev r q2
      helper (EQBranch ev q1    q2) r = EQBranch ev q2 (eqMerge q1 r)

eqInsert :: Event -> EventQueue -> EventQueue
eqInsert ev q = eqMerge q (eqSingle ev)

eqFromList :: [Event] -> EventQueue
eqFromList evs = foldr eqMerge EQNil . map eqSingle

eqRun :: (Event -> EventQueue -> ([a], EventQueue)) -> EventQueue -> [a]
eqRun _ EQNil = []
eqRun f q     = case uncurry f $ eqGet q of
                  (xs, q') -> xs ++ eqRun f q'

--- Sweep line
type SweepElemRef = IORef (Maybe SweepElem)
data SweepElem = SE !SegmentIndex !SweepElemRef !SweepElemRef
type SweepLine = (IM.IntMap Position,  SweepElem)

slEmpty :: SweepLine
slEmpty = IM.empty

slInsert :: SegmentIndex -> SweepLine -> IO (SweepLine, [Neighbors])
slInsert s sl = do
  newLeft  <- newIORef Nothing
  newRight <- newIORef Nothing
  let newElem = SE s newLeft newRight

  let update _    set | IM.null set = return []
      update left set = do
        let (extrm, oldElem@(SE _ l r)) = f set
                where f = if left then IM.findMax else IM.findMin
        let newRef = (if left then newLeft else newRight)
            oldRef = (if left then r else l)
        old <- readIORef oldRef
        writeIORef oldRef newElem
        writeIORef newRef oldElem
        return [extrm]

  let (ltSet, gtSet) = IM.split s sl
  l <- update True  ltSet
  r <- update False gtSet
  return (IM.insert s newElem sl, l ++ r)

slSwap :: Neighbors -> SweepLine -> (SweepLine, [Neighbors])
slSwap (s1,s2) sl =
  let (s1L, s1R) = sl IM.! s1
      (s2L, s2R) = sl IM.! s2

      newSl = IM.insert s1 (s2L, s2R) $
              IM.insert s2 (s1L, s1R) $ sl
      changes = (if s1L == siNone then id else ((s1L,s2):))
                (if s2R == siNone then [] else [(s1,s2R)])
  in if interDebug (s1R /= s2 || s2L /= s1)
     then interErr "[slSwap] they're not neighbors!"
     else (newSl, changes)




-- -- This looks like Data.Set, however with some peculiarities.
-- type SLSize = Int
-- data SweepLine = SLNil
--                | SLBranch !SLSize !SegmentIndex SweepLine SweepLine

-- slSize :: SweepLine -> SLSize
-- slSize SLNil                 = 0
-- slSize (SLBranch size _ _ _) = size

-- slSingle :: SegmentIndex -> SweepLine
-- slSingle s = SLBranch 1 s SLNil SLNil

-- slInsert :: SegmentIndex -> SweepLine -> SweepLine
-- slInsert s SLNil                 = slSingle s
-- slInsert s (SLBranch size t l r) =
--     case s `compare` t of
--       LT -> slBalance t (slInsert s l) r
--       GT -> slBalance t l (slInsert s r)
--       EQ -> interErr "[slInsert] segment already on sweepline"

-- slRemove :: SegmentIndex -> SweepLine -> SweepLine
-- slRemove s SLNil = interErr "[slRemove] segment not on sweepline"
-- slRemove s (SLBranch size t l r) =
--     case s `compare` t of
--       LT -> slBalance t (slRemove s l) r
--       GT -> slBalance t l (slRemove s r)

-- slBalance

--- Bentley-Ottmann functions
boEvents :: IndSeg -> [Event]
boEvents (index, a0,a1) = [EvStart start index, EvEnd end index]
    where (start,end) | a0 < a1   = (a0,a1)
                      | otherwise = (a1,a0)

bentleyOttmann :: [IndSeg] -> [InterIndexes]
bentleyOttmann isegs = eqRun go . eqFromList . concatMap boEvents
    where
      segArray = array (0, length isegs - 1) isegs

      go EQNil = []
      go


--------- END OF BENTLEY-OTTMANN ---------

-}