First, the new definitions:<\/p>\n
data Region = UnitCircle\r\n | Polygon [Coordinate]\r\n | AffineTransform Matrix3x3 Region\r\n | Empty\r\n deriving Show\r\n\r\ntype Vector3 = (Float, Float, Float)\r\ntype Matrix3x3 = (Vector3, Vector3, Vector3)\r\ntype Matrix2x2 = ((Float, Float), (Float, Float))<\/pre>\nAnd some old ones that will suffice unchanged:<\/p>\n
type Coordinate = (Float, Float)\r\ntype Ray = (Coordinate, Coordinate)\r\n\r\nisLeftOf :: Coordinate -> Ray -> Bool\r\n(px,py) `isLeftOf` ((ax,ay), (bx,by))\r\n = let (s,t) = (px-ax, py-ay)\r\n (u,v) = (px-bx, py-by)\r\n in s * v >= t * u\r\n\r\nisRightOf :: Coordinate -> Ray -> Bool\r\n(px,py) `isRightOf` ((ax,ay), (bx,by))\r\n = let (s,t) = (px-ax, py-ay)\r\n (u,v) = (px-bx, py-by)\r\n in s * v <= t * u\r\n\r\ncontainsR :: Region -> Coordinate -> Bool\r\n(Polygon pts) `containsR` p\r\n = let leftOfList = map isLeftOfp (zip pts (tail pts ++ [head pts]))\r\n isLeftOfp p' = isLeftOf p p'\r\n rightOfList = map isRightOfp (zip pts (tail pts ++ [head pts]))\r\n isRightOfp p' = isRightOf p p'\r\n in and leftOfList || and rightOfList\r\nEmpty `containsR` p = False<\/pre>\nNow we have the new forms of containsR<\/tt> to write. The UnitCircle<\/tt> is easy:<\/p>\nUnitCircle `containsR` (x,y) = x^2 + y^2 <= 1<\/pre>\nAnd in general, the transform of a region contains a point if the region contains the inverse transform of that point.<\/p>\n
(AffineTransform m r) `containsR` (x,y)\r\n = if determinant3 m == 0\r\n then singularContains m (x,y)\r\n else let m' = inverse m\r\n (x', y', _) = matrixMul m' (x,y,1)\r\n in r `containsR` (x', y')<\/pre>\nNow some standard code for multiplying and inverting matrices:<\/p>\n
matrixMul :: Matrix3x3 -> Vector3 -> Vector3\r\nmatrixMul (r1, r2, r3) v\r\n = (dotProduct r1 v,\r\n dotProduct r2 v,\r\n dotProduct r3 v)\r\n\r\ndotProduct :: Vector3 -> Vector3 -> Float\r\ndotProduct (a,b,c) (x,y,z) = a*x + b*y + c*z\r\n\r\ninverse :: Matrix3x3 -> Matrix3x3\r\ninverse ((a,b,c), (d,e,f), (g,h,i))\r\n = let det = determinant3 ((a,b,c), (d,e,f), (g,h,i))\r\n a' = determinant2 ((e,f), (h,i)) \/ det\r\n b' = determinant2 ((c,b), (i,h)) \/ det\r\n c' = determinant2 ((b,c), (e,f)) \/ det\r\n d' = determinant2 ((f,d), (i,g)) \/ det\r\n e' = determinant2 ((a,c), (g,i)) \/ det\r\n f' = determinant2 ((c,a), (f,d)) \/ det\r\n g' = determinant2 ((d,e), (g,h)) \/ det\r\n h' = determinant2 ((b,a), (h,g)) \/ det\r\n i' = determinant2 ((a,b), (d,e)) \/ det\r\n in ((a',b',c'), (d',e',f'), (g',h',i'))\r\n\r\ndeterminant3 :: Matrix3x3 -> Float\r\ndeterminant3 ((a,b,c), (d,e,f), (g,h,i))\r\n = let aei = a * e * i\r\n afh = a * f * h\r\n bdi = b * d * i\r\n cdh = c * d * h\r\n bfg = b * f * g\r\n ceg = c * e * g\r\n in aei - afh - bdi + cdh + bfg - ceg\r\n\r\ndeterminant2 :: Matrix2x2 -> Float\r\ndeterminant2 ((a,b), (c,d))\r\n = a * d - b * c<\/pre>\nThe only thing left is the question of how to deal with a singular (non-invertible) matrix. If an affine matrix is non-invertible, that means that AE - BD = 0. Either all four coefficients are 0, or AE = BD. If all 4 are 0, then every point in the region will be collapsed into the single point (C,F), and we need to check that this is the given point. If AE = BD otherwise, we have for any point (x,y)<\/tt> in the region:<\/p>\nx' = Ax + By + C\r\ny' = Dx + Ey + F\r\nDx' = ADx + BDy + CD\r\nAy' = ADx + AEy + AF\r\n = ADx + BDy + AF [AE = BD]\r\n\u00e2\u02c6\u00b4 Dx' - Ay' = CD - AF\r\nDx' - Ay' + AF - CD = 0<\/pre>\nso any point in the region will be collapsed into the line Dx' - Ay' + AF - CD = 0, and we need to check that the given point satisfies this equation.<\/p>\n
singularContains :: Matrix3x3 -> Coordinate -> Bool\r\nsingularContains ((a,b,c), (d,e,f), _) (x,y)\r\n = if (a == 0 && b == 0 && d == 0 && e == 0)\r\n then (x == c && y == f)\r\n else (d*x - a*y + a*f - c*d == 0)<\/pre>\n","protected":false},"excerpt":{"rendered":"First, the new definitions: data Region = UnitCircle | Polygon [Coordinate] | AffineTransform Matrix3x3 Region | Empty deriving Show type Vector3 = (Float, Float, Float) type Matrix3x3 = (Vector3, Vector3, Vector3) type Matrix2x2 = ((Float, Float), (Float, Float)) And some old ones that will suffice unchanged: type Coordinate = (Float, Float) type Ray = (Coordinate, […]<\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"open","ping_status":"open","sticky":false,"template":"","format":"standard","meta":[],"categories":[1],"tags":[],"_links":{"self":[{"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=\/wp\/v2\/posts\/42"}],"collection":[{"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=%2Fwp%2Fv2%2Fcomments&post=42"}],"version-history":[{"count":0,"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=\/wp\/v2\/posts\/42\/revisions"}],"wp:attachment":[{"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=%2Fwp%2Fv2%2Fmedia&parent=42"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=%2Fwp%2Fv2%2Fcategories&post=42"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/www.elbeno.com\/haskell_soe_blog\/index.php?rest_route=%2Fwp%2Fv2%2Ftags&post=42"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}