Difference between revisions of "Graham Scan Implementation"
Jump to navigation
Jump to search
(Add test code) |
(make the implementation neat) |
||
Line 5: | Line 5: | ||
import Data.Function (on) |
import Data.Function (on) |
||
− | + | data Direction = LEFT | RIGHT | STRAIGHT |
|
⚫ | |||
− | data Direction = GoLeft | GoRight | GoStraight |
||
+ | |||
⚫ | |||
⚫ | |||
− | -- determine direction change via point a->b->c |
||
− | -- -- define 2D point |
||
⚫ | |||
isTurned :: Pt -> Pt -> Pt -> Direction |
isTurned :: Pt -> Pt -> Pt -> Direction |
||
isTurned (Pt (ax, ay)) (Pt (bx, by)) (Pt (cx, cy)) = case sign of |
isTurned (Pt (ax, ay)) (Pt (bx, by)) (Pt (cx, cy)) = case sign of |
||
− | + | EQ -> STRAIGHT |
|
− | + | LT -> RIGHT |
|
− | + | GT -> LEFT |
|
− | + | where sign = compare ((bx - ax) * (cy - ay)) ((cx - ax) * (by - ay)) |
|
⚫ | |||
− | |||
− | -- implement Graham scan algorithm |
||
− | |||
− | -- -- Helper functions |
||
⚫ | |||
⚫ | |||
− | buttonLeft [] = Pt (1/0, 1/0) |
||
− | buttonLeft [pt] = pt |
||
− | buttonLeft (pt:pts) = minY pt (buttonLeft pts) where |
||
− | minY (Pt (ax, ay)) (Pt (bx, by)) |
||
⚫ | |||
− | | ay < by = Pt (ax, ay) |
||
− | | ax < bx = Pt (ax, ay) |
||
− | | otherwise = Pt (bx, by) |
||
− | |||
− | -- -- Main |
||
− | convex :: [Pt] -> [Pt] |
||
− | convex [] = [] |
||
− | convex [pt] = [pt] |
||
− | convex [pt0, pt1] = [pt0, pt1] |
||
− | convex pts = scan [pt0] spts where |
||
− | -- Find the most buttonleft point pt0 |
||
− | pt0 = buttonLeft pts |
||
⚫ | |||
− | -- Sort other points `ptx` based on angle <pt0->ptx> |
||
+ | gScan pts |
||
⚫ | |||
+ | | length pts >= 3 = scan [pt0] rests |
||
⚫ | |||
+ | | otherwise = pts |
||
− | {-the secondary key make sure collinear points in order-} |
||
+ | where |
||
⚫ | |||
⚫ | |||
+ | pt0 = foldr bottomLeft (Pt (1/0, 1/0)) pts where |
||
+ | bottomLeft pa pb = case ord of |
||
⚫ | |||
⚫ | |||
⚫ | |||
⚫ | |||
− | -- |
+ | -- Sort other points based on angle |
⚫ | |||
⚫ | |||
⚫ | |||
⚫ | |||
− | + | abs (x - x0)) |
|
− | where pz = last ps |
||
+ | -- Scan the points to find out convex |
||
⚫ | |||
⚫ | |||
⚫ | |||
⚫ | |||
⚫ | |||
− | + | | isTurned pz p0 p1 == STRAIGHT = [pz, p0] |
|
− | + | where pz = last ps |
|
− | |||
− | |||
− | -- Test data |
||
− | pts1 = [Pt (0,0), Pt (1,0), Pt (2,1), Pt (3,1), Pt (2.5,1), Pt (2.5,0), Pt (2,0)] |
||
⚫ | |||
⚫ | |||
− | -- ---- v |
||
− | -- (0,0)------->(1,0)-- (2,0)<-(2.5,0) |
||
− | pts2 = [Pt (1,-2), Pt (-1, 2), Pt (-3, 6), Pt (-7, 3)] |
||
− | -- -- This case demonstrate a collinear points across the button-left point |
||
− | -- |-----(-3,6) |
||
− | -- |----- ^__ |
||
− | -- <----- |--- |
||
− | -- (-7,3) | |
||
− | -- |-- (-1,2) |
||
− | -- |---- <-- |
||
− | -- |------------- |--| |
||
− | -- |------------>(1,-2) |
||
− | pts3 = [Pt (6,0), Pt (5.5,2), Pt (5,2), Pt (0,4), Pt (1,4), Pt (6,4), Pt (0,0)] |
||
− | -- -- This case demonstrates a consecutive inner points |
||
− | -- (0,4)<-(1,4)<--------------------(6,4) |
||
− | -- | |-----------^ |
||
⚫ | |||
− | -- (0,0)--------------------------->(6,0) |
||
⚫ | |||
− | pts4 = [Pt (0,0), Pt (0,4), Pt (1,4), Pt (2,3), Pt (6,5), Pt (6,0)] |
||
⚫ | |||
− | -- -- This case demonstrates a consecutive inner points |
||
⚫ | |||
⚫ | |||
⚫ | |||
− | -- collinear points |
||
− | pts5 = [Pt (5,5), Pt (4,4), Pt (3,3), Pt (2,2), Pt (1,1)] |
||
− | pts6 = reverse pts5 |
||
</haskell> |
</haskell> |
||
Revision as of 14:30, 5 January 2015
Descriptions of this problem can be found in Real World Haskell, Chapter 3
import Data.List (sortBy)
import Data.Function (on)
data Direction = LEFT | RIGHT | STRAIGHT
deriving (Show, Eq)
data Pt = Pt (Double, Double)
isTurned :: Pt -> Pt -> Pt -> Direction
isTurned (Pt (ax, ay)) (Pt (bx, by)) (Pt (cx, cy)) = case sign of
EQ -> STRAIGHT
LT -> RIGHT
GT -> LEFT
where sign = compare ((bx - ax) * (cy - ay)) ((cx - ax) * (by - ay))
gScan :: [Pt] -> [Pt]
gScan pts
| length pts >= 3 = scan [pt0] rests
| otherwise = pts
where
-- Find the most bottom-left point pt0
pt0 = foldr bottomLeft (Pt (1/0, 1/0)) pts where
bottomLeft pa pb = case ord of
LT -> pa
GT -> pb
EQ -> pa
where ord = (compare `on` (\ (Pt (x, y)) -> (y, x))) pa pb
-- Sort other points based on angle
rests = tail (sortBy (compare `on` compkey pt0) pts) where
compkey (Pt (x0, y0)) (Pt (x, y)) = (atan2 (y - y0) (x - x0),
abs (x - x0))
-- Scan the points to find out convex
-- -- handle the case that all points are collinear
scan [p0] (p1:ps)
| isTurned pz p0 p1 == STRAIGHT = [pz, p0]
where pz = last ps
scan (x:xs) (y:z:rsts) = case isTurned x y z of
RIGHT -> scan xs (x:z:rsts)
STRAIGHT -> scan (x:xs) (z:rsts) -- skip collinear points
LEFT -> scan (y:x:xs) (z:rsts)
scan xs [z] = z : xs
Test
import Test.QuickCheck
import Control.Monad
import Data.List
convex' = map (\(Pt x) -> x) . convex . map Pt
prop_convex n = forAll (replicateM n arbitrary) (\xs -> sort (convex' xs) == sort (convex' (convex' xs)))
The results are:
> quickCheck (prop_convex 0)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 1)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 2)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 3)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 4)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 5)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 10)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 100)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 1000)
+++ OK, passed 100 tests.
> quickCheck (prop_convex 10000)
+++ OK, passed 100 tests.