path: root/src/Geometry.hs

module Geometry
  ( Line (..),
    lineThrough2Points,
    lineThroughPointWithAngle,
    LineSide (..),
    sideOfLine,
    proj,
    perp',
    orthoproj,
    distPointLine,
    distPointSegment,
    rotate,
    -- Rectangle types
    GeometryRectangle (..),
    rectangleToGeometry,
    geometryToRectangle,
    -- Intersection
    IntersectionResult (..),
    getIntersectingRectangles,
  )
where

import Data.Int (Int32)
import Data.List (sortOn)
import Data.Maybe (catMaybes, listToMaybe)
import Graphics.X11 (Rectangle (..))
import Prelude hiding (subtract)

-- | A line in parametric form: (x + t*dx, y + t*dy) for t ∈ ℝ
data Line a = Line
  { x :: !a,
    y :: !a,
    dx :: !a,
    dy :: !a
  }
  deriving (Eq, Show, Read)

-- | Construct a line passing through two points
lineThrough2Points :: (Num a, Ord a) => (a, a) -> (a, a) -> Line a
lineThrough2Points (x1, y1) (x2, y2) =
  Line
    { x = x1,
      y = y1,
      dx = x2 - x1,
      dy = y2 - y1
    }

-- | Construct a line through a point with a given angle (in radians)
--   Angle is measured counterclockwise from the positive x-axis
lineThroughPointWithAngle :: (Floating a, Ord a) => (a, a) -> a -> Line a
lineThroughPointWithAngle (x, y) angle =
  Line
    { x = x,
      y = y,
      dx = cos angle,
      dy = sin angle
    }

-- | Which side of the line is a point on?
data LineSide
  = LeftSide
  | RightSide
  | OnLine
  deriving (Eq, Show)

-- | Determine which side of the line a point lies on
--   Uses the cross product to determine orientation
sideOfLine :: (Fractional a, Ord a) => Line a -> (a, a) -> LineSide
sideOfLine line (px, py) =
  case compare cross 0 of
    GT -> LeftSide
    LT -> RightSide
    EQ -> OnLine
  where
    -- Vector from line origin to point
    relX = px - x line
    relY = py - y line
    -- Cross product in 2D: (dx, dy) × (px, py) = dx*py - dy*px
    cross = dx line * relY - dy line * relX

-- | Project point onto line
proj :: (Fractional a, Ord a) => Line a -> (a, a) -> (a, a)
proj line (px, py) =
  (x line + t * dx line, y line + t * dy line)
  where
    -- Vector from line origin to point
    relX = px - x line
    relY = py - y line
    -- t = (rel · direction) / (direction · direction)
    dotProd = relX * dx line + relY * dy line
    dirSq = dx line * dx line + dy line * dy line
    t = dotProd / dirSq

-- | Perpendicular vector from point to line
perp' :: (Fractional a, Ord a) => Line a -> (a, a) -> (a, a)
perp' line point =
  (px - projX, py - projY)
  where
    (projX, projY) = proj line point
    (px, py) = point

-- | Orthogonal projection (same as proj, but returns the projected point)
orthoproj :: (Fractional a, Ord a) => Line a -> (a, a) -> (a, a)
orthoproj = proj

-- | Distance from point to line
distPointLine :: (Floating a, Ord a) => Line a -> (a, a) -> a
distPointLine line point =
  sqrt (dx * dx + dy * dy)
  where
    (dx, dy) = perp' line point

-- | Distance from point to line segment (not infinite line)
distPointSegment :: (Floating a, Ord a) => (a, a) -> (a, a) -> (a, a) -> a
distPointSegment (ax, ay) (bx, by) (px, py)
  | t <= 0 = sqrt (dx1 * dx1 + dy1 * dy1)
  | t >= 1 = sqrt (dx2 * dx2 + dy2 * dy2)
  | otherwise = sqrt (dx * dx + dy * dy)
  where
    abX = bx - ax
    abY = by - ay
    apX = px - ax
    apY = py - ay
    t = (apX * abX + apY * abY) / (abX * abX + abY * abY)
    projX = ax + t * abX
    projY = ay + t * abY
    dx = px - projX
    dy = py - projY
    dx1 = px - ax
    dy1 = py - ay
    dx2 = px - bx
    dy2 = py - by

-- | Rotate a point around the origin by an angle (in radians)
rotate :: (Floating a, Ord a) => a -> (a, a) -> (a, a)
rotate angle (x, y) =
  (x * cos angle - y * sin angle, x * sin angle + y * cos angle)

-- | A simple Rectangle type for geometry operations
--   Uses (x, y) as the top-left corner, with width and height
data GeometryRectangle a = GeometryRectangle
  { geoX :: !a,
    geoY :: !a,
    geoWidth :: !a,
    geoHeight :: !a
  }
  deriving (Eq, Show, Read)

-- | Convert X11's Rectangle to Geometry's Rectangle
rectangleToGeometry :: Num a => Graphics.X11.Rectangle -> GeometryRectangle a
rectangleToGeometry (Graphics.X11.Rectangle x y w h) =
  GeometryRectangle
    { geoX = fromIntegral x,
      geoY = fromIntegral y,
      geoWidth = fromIntegral w,
      geoHeight = fromIntegral h
    }

-- | Convert Geometry's Rectangle to X11's Rectangle
geometryToRectangle :: Integral a => GeometryRectangle a -> Graphics.X11.Rectangle
geometryToRectangle (GeometryRectangle x y w h) =
  Graphics.X11.Rectangle (fromIntegral x) (fromIntegral y) (fromIntegral w) (fromIntegral h)

-- | Alias for GeometryRectangle for backwards compatibility
type Rectangle a = GeometryRectangle a

-- | Result of intersecting a line with a rectangle
data IntersectionResult a b = IntersectionResult
  { tag :: a,
    rectangle :: GeometryRectangle b,
    intersectionPoint :: (b, b)
  }
  deriving (Eq, Show)

-- | Find all rectangles that intersect the given line, sorted by distance from origin
--   Returns tags of intersected rectangles in order of their intersection point
--   closest to (0,0)
getIntersectingRectangles :: (Floating a, Ord a) => Line a -> [(GeometryRectangle a, b)] -> [IntersectionResult b a]
getIntersectingRectangles line rects = catMaybes $ map (checkIntersection line) rects

-- | Check if a line intersects a rectangle, and if so, return the intersection point
--   closest to the origin
checkIntersection :: (Floating a, Ord a) => Line a -> (GeometryRectangle a, b) -> Maybe (IntersectionResult b a)
checkIntersection line (rect, tag) = IntersectionResult tag rect <$> findClosestIntersection line rect

-- | Find the intersection point closest to the origin
findClosestIntersection :: (Floating a, Ord a) => Line a -> GeometryRectangle a -> Maybe (a, a)
findClosestIntersection line rect =
  listToMaybe $ sortOn distanceFromOrigin $ filter (pointInRectangle rect) candidates
  where
    -- Find intersections with each edge of the rectangle
    topEdge = lineThrough2Points (geoX rect, geoY rect) (geoX rect + geoWidth rect, geoY rect)
    bottomEdge = lineThrough2Points (geoX rect, geoY rect + geoHeight rect) (geoX rect + geoWidth rect, geoY rect + geoHeight rect)
    leftEdge = lineThrough2Points (geoX rect, geoY rect) (geoX rect, geoY rect + geoHeight rect)
    rightEdge = lineThrough2Points (geoX rect + geoWidth rect, geoY rect) (geoX rect + geoWidth rect, geoY rect + geoHeight rect)

    candidates =
      catMaybes
        [ intersectLines line topEdge,
          intersectLines line bottomEdge,
          intersectLines line leftEdge,
          intersectLines line rightEdge
        ]

-- | Check if a point is inside a rectangle (including edges)
pointInRectangle :: (Num a, Ord a) => GeometryRectangle a -> (a, a) -> Bool
pointInRectangle rect (px, py) =
  px >= geoX rect
    && px <= geoX rect + geoWidth rect
    && py >= geoY rect
    && py <= geoY rect + geoHeight rect

-- | Find the intersection point of two lines, if it exists
intersectLines :: (Fractional a, Eq a) => Line a -> Line a -> Maybe (a, a)
intersectLines (Line x1 y1 dx1 dy1) (Line x2 y2 dx2 dy2) =
  -- Solve: x1 + t*dx1 = x2 + s*dx2 and y1 + t*dy1 = y2 + s*dy2
  -- Using Cramer's rule
  let det = dx1 * (- dy2) - dy1 * (- dx2)
      detT = (x2 - x1) * (- dy2) - (y2 - y1) * (- dx2)
   in if det == 0
        then Nothing -- Parallel lines
        else
          let t = detT / det
              intersectionX = x1 + t * dx1
              intersectionY = y1 + t * dy1
           in Just (intersectionX, intersectionY)

-- | Euclidean distance from origin
distanceFromOrigin :: (Num a, Floating a) => (a, a) -> a
distanceFromOrigin (x, y) = sqrt (x * x + y * y)