photo filter now uses voroni diagram

There's still a couple bugs, some of the shapes are kinda long and
pointy. In the meantime, shapes are placed on the board longest
perimeter first to cover it up.

Also, there's still some gaps occasionally.

image-triangles.cabal

@@ -82,6 +82,7 @@ executable image-triangles
                      , diagrams-lib
                      , diagrams-cairo
                      , diagrams-svg
+                     , diagrams-contrib
                      , parallel
                      , linear
                      , vector
@@ -102,4 +103,3 @@ executable image-triangles
                "-with-rtsopts= -N"
               -- -prof
 --               -fexternal-interpreter
---               import Servant.Client

src/Main.hs

@@ -32,6 +32,9 @@ import Control.Arrow
 import Data.Colour.RGBSpace (uncurryRGB)
 import qualified Control.Monad.Parallel as MP
 import qualified CircumCircle as CC
+import qualified Debug.Trace as T
+import qualified Diagrams as DP
+import qualified Diagrams.Prelude as D
 
 data Options = Options {
                     numPoints :: Int,
@@ -53,52 +56,70 @@ convImage = Vec.map tosRGB' . Int.toVector
 corners :: [(Double, Double)]
 corners = (,) <$> [0, 1] <*> [0, 1]
 
-scalePointToImage :: (Int, Int) -> Point V2 Double -> Point V2 Int
+-- scalePointToImage :: (Int, Int) -> Point V2 Double -> Point V2 Int
-scalePointToImage (ymax, xmax) p = round <$> (p2 (fromIntegral xmax, fromIntegral ymax) * p)
+-- scalePointToImage (ymax, xmax) p = round <$> (p2 (fromIntegral xmax, fromIntegral ymax) * p)
 
-genImage :: Image VU G.RGB Double -> StdGen -> Int -> QDiagram B V2 Double Any
+-- genImage :: Image VU G.RGB Double -> StdGen -> Int -> QDiagram B V2 Double Any
-genImage image gen cornerCount =
+-- genImage image gen cornerCount =
-  scaleY widthHeightRatio
+--   scaleY widthHeightRatio
-  . reflectY
+--   . reflectY
-  . rectEnvelope (mkP2 0 0) (1 ^& 1)
+--   . rectEnvelope (mkP2 0 0) (1 ^& 1)
-  . mconcat
+--   . mconcat
-  . map (uncurry Ren.placeTri . second (uncurryRGB CL.rgb))
+--   . map (uncurry Ren.placeTri . second (uncurryRGB CL.rgb))
-  . withStrategy (parListChunk 1000 rdeepseq)
+--   . withStrategy (parListChunk 1000 rdeepseq)
-  $ My.mapMaybe (\tri ->  (tri,) . CL.toRGB <$> triColor tri) triangles
+--   $ My.mapMaybe (\tri ->  (tri,) . CL.toRGB <$> triColor tri) triangles
-  where
+--   where
-    widthHeightRatio = (fromIntegral $ fst dimensions) / ((fromIntegral $ snd dimensions) :: Double)
+--     widthHeightRatio = (fromIntegral $ fst dimensions) / ((fromIntegral $ snd dimensions) :: Double)
 
-    img' = convImage image
+--     img' = convImage image
-    dimensions = Img.dims image
+--     dimensions = Img.dims image
-    triangles = S.toList . Tri.findTriangles . withStrategy (parListChunk 1000 rdeepseq) . Tri.toPlanarGraph . take cornerCount . map p2 $ corners ++ Tri.randomPoints gen
+--     triangles = S.toList . Tri.findTriangles . withStrategy (parListChunk 1000 rdeepseq) . Tri.toPlanarGraph . take cornerCount . map p2 $ corners ++ Tri.randomPoints gen
 
-    triColor tri = getTriangleAverageRGB img' dimensions <$> (if S.size scaled == 3 then Just scaled else Nothing)
+--     triColor tri = getTriangleAverageRGB img' dimensions <$> (if S.size scaled == 3 then Just scaled else Nothing)
-      where
+--       where
-        scaled = S.map (scalePointToImage dimensions) tri
+--         scaled = S.map (scalePointToImage dimensions) tri
 
 
+shapeCircumference = Data.List.sum . map D.norm . loopOffsets . fromVertices
 
--- genImage' :: Image VU G.RGB Double -> StdGen -> Int -> QDiagram Cairo V2 Double Any
+genImage' :: Image VU G.RGB Double -> StdGen -> Int -> QDiagram SVG V2 Double Any
 genImage' image gen cornerCount =  -- image gen cornerCount =
   scaleY widthHeightRatio
   . reflectY
   . rectEnvelope (mkP2 0 0) (1 ^& 1)
-  . position
+  -- . atop visualizeGraph
-  . map (\shape -> (head shape, ) . showOrigin
+  . mconcat
-                    . fillColor (Tri.getShapeAverageRGB img' dimensions . map (scalePointToImage dimensions) $ shape)
+  . map drawVoroniRegion
-                    . strokeLoop  . closeLine . fromVertices $ shape)
+  -- . concatMap (\(center, shape) -> (zip shape (repeat $ circle 0.002)))
-  -- . withStrategy (parListChunk 1000 rdeepseq)
+  . sortOn shapeCircumference
+  . withStrategy (parListChunk 50 rdeepseq)
+  . map (uncurry Tri.voroniDiagramCorners)
   $ voroni
   where
-    widthHeightRatio = (fromIntegral $ fst dimensions) / ((fromIntegral $ snd dimensions) :: Double)
 
-    triColor tri = getTriangleAverageRGB img' dimensions <$> (if S.size scaled == 3 then Just scaled else Nothing)
+    drawVoroniRegion shape = lw 0
-      where
+                    . fillColor (Tri.voroniRegionAverageColor img' dimensions shape )
-        scaled = S.map (scalePointToImage dimensions) tri
+                    . strokeLocLoop . fromVertices $ shape
+
+
+    widthHeightRatio = (fromIntegral $ fst dimensions) / ((fromIntegral $ snd dimensions) :: Double)
 
     img' = convImage image
     dimensions = Img.dims image
-    voroni = take 40 . Tri.findVoroniDiagram . withStrategy (parListChunk 1000 rdeepseq) . Tri.toPlanarGraph . take cornerCount . map p2 $ corners ++ Tri.randomPoints gen
+
+    singleVoroni = last voroni
+
+    visualizeGraph :: QDiagram SVG V2 Double Any
+    visualizeGraph = lc red . lw 1 . position
+      . map (\(x0, x1) -> (x0,) . strokeLine $ (x0 ~~ x1))
+      . Tri.toPlanarGraph
+      $ corners'
+
+    voroni = -- take 14 . drop 24 .
+      Tri.findVoroniDiagram -- . withStrategy (parListChunk 1000 rdeepseq)
+      . Tri.toPlanarGraph $ corners'
+
+    corners' = take cornerCount . map p2 $ corners ++ Tri.randomPoints gen
 
 deriving instance Generic (CL.RGB a)
 deriving instance NFData a => NFData (CL.RGB a)
@@ -118,11 +139,11 @@ main :: IO ()
 main = do
   CLIOptions{..} <- getRecord "image options"
   let (Options {gen = gen}) = defaultOpts
-  gen' <- pure . mkStdGen $ 2345 -- getStdGen
+  gen' <- getStdGen -- pure . mkStdGen $ 2344
   let gens :: [StdGen] = map fst . iterate (split . snd) . split $ gen'
   print gen'
   image <- Img.readImageRGB VU input
-  let diagram = genImage image gen' cornerCount
+  -- let diagram = genImage image gen' cornerCount
   let nums = zip gens $ map show [0..60]
   let dimVector = toDimensionVector image
   renderSVG output dimVector (genImage' image gen' cornerCount)

src/Triangles.hs

@@ -29,6 +29,14 @@ import qualified Data.Map as M
 import qualified Data.Set as S
 import Control.Arrow
 import qualified Data.List as L
+import qualified Data.Function as F
+import qualified Diagrams.TwoD.Path.IntersectionExtras as I
+import Diagrams.Trail (trailPoints)
+import qualified Graphics.Image as H
+import Debug.Trace (traceShow)
+import qualified Data.Colour.Names as C
+import qualified Debug.Trace as D
+import Diagrams.TwoD.Segment.Bernstein (listToBernstein)
 -- import qualified Linear.Affine as L
 type Image_ = Vec.Vector Pixel_
 type Pixel_ = Colour Double
@@ -53,45 +61,110 @@ randomPoints :: StdGen -> [(Double, Double)]
 randomPoints = map (bimap toZeroToOneTuple toZeroToOneTuple) . randomRs ((0 :: Word, 0 :: Word), (maxBound, maxBound))
   where
     toZeroToOneTuple :: Word -> Double
-    toZeroToOneTuple x = ((fromIntegral x / (fromIntegral (maxBound :: Word))) * 1.1) - 0.05
+    toZeroToOneTuple x = ((fromIntegral x / (fromIntegral (maxBound :: Word))) * 1.2) - 0.1
 
 combinations :: (Ord b, Floating b, NFData b) => [P2 b] -> [(P2 b, P2 b)]
 combinations xs =
   sortOn (abs . uncurry distanceA)
   . S.toList . S.fromList
   . filter (uncurry (/=))
-  . concat . withStrategy (parListChunk 2000 rdeepseq)
+  . sortOn (abs . uncurry distanceA)
-  . map (\x -> take edgeLengthThreshold . sortOn (abs . uncurry distanceA)
+  -- . concat . withStrategy (parListChunk 2000 rdeepseq) $
-  . map (\y -> (min x y, max x y)) $ xs) $ xs
+   $ [(x,y) | x <- xs, y <- xs]
+  -- . map (\x -> take edgeLengthThreshold .
+  --         sortOn (abs . uncurry distanceA)
+  -- . map (\y -> (min x y, max x y)) $ xs) $ xs
   where
-    edgeLengthThreshold = 45
+    edgeLengthThreshold = 10
 
--- toPlanarGraph :: (NFData n, Floating n, Ord n) => [Point V2 n] -> [(Point V2 n, Point V2 n)]
+toPlanarGraph :: forall n . (NFData n, Floating n, Ord n) => [P2 n] -> [(Point V2 n, Point V2 n)]
-toPlanarGraph :: (NFData n, Floating n, Ord n) => [P2 n] -> [(Point V2 n, Point V2 n)]
 toPlanarGraph points =
-  removeIntersections . combinations $ points
+  removeIntersections .
+  combinations $ points
   where
     numPoints = length points
 
+    removeIntersections :: [(Point V2 n, Point V2 n)] -> [(Point V2 n, Point V2 n)]
     removeIntersections = foldl' addIfNoIntersection []
 
-    addIfNoIntersection xs x = if all (noIntersection x) xs then (x:xs) else xs
 
-    noIntersection l1 l2 = (==) sharedEndPoint . length $ intersectPointsT (uncurry toLocatedTrail l1) (uncurry toLocatedTrail l2)
+    addIfNoIntersection xs x | all (noIntersection x) xs =  x:xs
+                             | otherwise = xs
+
+    noIntersection l1 l2 = sharedEndPoint || (null $ intersectPointsT (uncurry toLocatedTrail l1) (uncurry toLocatedTrail l2))
       where
-        sharedEndPoint = (-) 4 . length . nub $ [fst l1, snd l1, fst l2, snd l2]
+        sharedEndPoint = (< 4) . length . nub $ [fst l1, snd l1, fst l2, snd l2]
 
 toLocatedTrail :: TrailLike a => Point (V a) (N a) -> Point (V a) (N a) -> Located a
 toLocatedTrail p1 p2 = fromVertices [p1, p2] `at` p1
 
+-- data VoroniRegion = VoroniRegion {
+--   center :: P2 Double,
+--   neighbors :: [P2 Double]
+--                                  }
+
+-- voroniDiagramCorners :: Point V2 Double -> [Point V2 Double] -> [Point V2 Double]
+-- voroniDiagramCorners :: Point (V c) (N c) -> [Point V2 Double] -> c
+-- voroniDiagramCorners :: forall t. (N t ~ Double, V t ~ V2, TrailLike t) => Point V2 Double -> [Point V2 Double] -> t
+withinShape :: RealFloat v => Point V2 v -> [Point V2 v] -> Point V2 v -> Bool
+withinShape pointInShape verticies candidate = all ((< quarterTurn) . fmap abs . uncurry signedAngleBetweenDirs) $ zip (shapeDirections pointInShape) (shapeDirections candidate)
+  where
+    shapeDirections p = map (dirBetween p) verticies
+
+
+sortOnAngle center = sortOn (normalizeAngle . signedAngleBetweenDirs xDir . dirBetween center)
+
+-- voroniDiagramCorners :: forall n . (Floating n, Real n, Show n) => Point V2 n -> [P2 n] -> [P2 n]
+voroniDiagramCorners center midpoints -- = midpoints
+  = sortOnAngle center . filter isValidMidpoint . concat $ [intersectPointsT l0 l1 | l0 <- tangentTrails, l1 <- tangentTrails]
+  where
+
+    lessThanQuarterTurn candidate = candidate <= (10001 / 40000) @@ turn || candidate >= (29999 / 40000) @@ turn
+      -- where
+      --   candidate' = (candidate ^. turn)  - (fromIntegral . floor $ (candidate ^. turn))
+    -- candidateVertecies = sortOn (normalizeAngle . signedAngleBetweenDirs xDir . dirBetween center) . concat $ [ intersectPointsT (tangentTrail x) (tangentTrail y)  | x <- midpoints, y <- midpoints ]
+
+    tangentTrails = map tangentTrail midpoints
+    -- validMidpoints = sortOn (normalizeAngle . signedAngleBetweenDirs xDir . dirBetween center) . filter isValidMidpoint $ midpoints
+
+    appendHead (x:xs) = xs ++ [x]
+
+    isValidMidpoint candidate = all isNonObtuseMidpoint . filter (/= candidate) $ midpoints
+      where
+        isNonObtuseMidpoint m = -- traceShow ((dirBetween m center, dirBetween m candidate), (candidate, m)) .
+          lessThanQuarterTurn . normalizeAngle $ angleBetweenDirs (dirBetween m center) (dirBetween m candidate)
+
+    tangentTrail midpoint = fromVertices $ [midpoint .-^ tangentVec, midpoint .+^ tangentVec]
+      where
+        -- implicitly uses the unit vector as an infinitely long vector
+        tangentVec = scale 8 .
+          fromDirection . rotateBy (1 / 4) $ dirBetween midpoint center
+
+
 -- findVoroniDiagram :: (Ord (v n), Additive v, Fractional n) => [(Point v n, Point v n)] -> [[Point v n]]
-findVoroniDiagram :: (Ord n, Ord (v n), Metric v, Floating n, R1 v, Real n) => [(Point v n, Point v n)] -> [[Point v n]]
+-- findVoroniDiagram :: (Ord n, Ord (v n), Metric v, Floating n, R1 v, Real n, r) => [(Point v n, Point v n)] -> [[Point v n]]
-findVoroniDiagram edges = M.elems . M.mapWithKey (\key -> L.sortOn (normalizeAngle . angleBetweenDirs xDir . dirBetween key) .
+findVoroniDiagram :: RealFloat n => [(Point V2 n, Point V2 n)] -> [(Point V2 n, [Point V2 n])]
-                                                   map (pointBetween key) . S.toList) $ adjacencyMap
+findVoroniDiagram edges =
+  M.toList
+  . M.mapWithKey (\key ->  L.sortOn (normalizeAngle . signedAngleBetweenDirs xDir . dirBetween key) . map (pointBetween key) .
+                   S.toList)
+  $ adjacencyMap
   where
     adjacencyMap = adjacencyMapOf edges
 
-    pointBetween p0 p1 = p0 .+^ ((p0 .-. p1) ^/ 2)
+    -- I'm not sure this part works
+    pointBetween p0 p1 = p0 .+^ ((p1 .-. p0) ^/ 2)
+
+-- filterNotInVoroniRegion :: RealFloat n => Point V2 n -> [Point V2 n] -> [Point V2 n] -> [Point V2 n]
+-- filterNotInVoroniRegion center midpoints = id -- filter allowed
+--   where
+--     allowed point
+--       = all ((< 0.25) . abs . (^. turn))
+--       . zipWith signedAngleBetweenDirs midpointAngles
+--       . map (dirBetween point) $ midpoints
+
+--     midpointAngles = map (dirBetween center) midpoints
+
 
 findTriangles :: Ord b => [(b, b)] -> S.Set (S.Set b)
 findTriangles edges = S.unions . S.map threeCyclesOf . M.keysSet $ adjacencyMap
@@ -119,15 +192,57 @@ getPointsInTriangle image pts = S.toList . S.unions . map S.fromList $ [
     where
         [p1, p2, p3] = sortOn fst . map (\(y, x) -> (x,y)) . map unp2 . S.toList $ pts
 
-blendEqually colors = C.affineCombo (zip (repeat fraction) colors) $ C.black
+blendEqually colors = C.affineCombo (zip (repeat fraction) colors) $ C.blue
     where
         fraction = 1.0 / (fromIntegral . length $ colors)
 
-getShapeAverageRGB :: Image_ -> (Int, Int) -> [P2 Int] -> C.Colour Double
+-- voroniRegionAverageColor :: Image_ -> (Int, Int) -> (P2 Int, [P2 Int]) -> C.Colour Double
-getShapeAverageRGB image sizes = blendEqually . concatMap (getColorsInTriangle image sizes) . filter ((== 3) . S.size) . map (S.fromList . take 3) . tails . L.sortOn (fst . unp2)
+voroniRegionAverageColor image (x', y') verticies
+  = blendEqually
+  . concatMap (getColorsInTriangle image (x', y'))
+  . filter ((== 3) . S.size)
+  . map (S.fromList . take 3)
+  . tails
+  . L.nub
+  . map scaleToImageCoords
+  $ verticies
+
+  where
+
+    scaleToImageCoords :: P2 Double -> P2 Int
+    scaleToImageCoords p = fmap round $ p2 ((fromIntegral x' * p ^. _x), fromIntegral y' * p ^. _y)
+
+    scaleToUnitCoords :: P2 Int -> P2 Double
+    scaleToUnitCoords p = p2 ((fromIntegral x' / (fromIntegral $ p ^. _x)), fromIntegral y' / (fromIntegral $ p ^. _y))
+  -- = blendEqually
+  -- . mapMaybe (index' . unp2)
+  -- . filter withinShape'
+  -- $ candidatePoints
+  -- where
+
+  --   withinShape' :: P2 Int -> Bool
+  --   withinShape' = withinShape center verticies . scaleToUnitCoords
+
+  --   candidatePoints = [p2 (x, y) | x <- [minX .. maxX], y <- [minY .. maxY]]
+
+  --   maxX = fst . unp2 $ maximumBy (compare `F.on` fst . unp2) verticies'
+  --   minX = fst . unp2 $ minimumBy (compare `F.on` fst . unp2) verticies'
+  --   maxY = snd . unp2 $ maximumBy (compare `F.on` snd . unp2) verticies'
+  --   minY = snd . unp2 $ minimumBy (compare `F.on` snd . unp2) verticies'
+
+  --   verticies' = map scaleToImageCoords verticies
+
+  --   index' :: (Int, Int) -> Maybe Pixel_
+  --   index' (y, x)
+  --       | y >= y' = Nothing
+  --       | x >= x' = Nothing
+  --       | y < 0   = Nothing
+  --       | x < 0   = Nothing
+  --       | otherwise = image Vec.!? ((y * x') + x)
+
 
 getColorsInTriangle :: Image_ -> (Int, Int) -> S.Set (P2 Int)  -> [C.Colour Double]
-getColorsInTriangle image (y', x') triangle = pixels
+getColorsInTriangle image (x', y') triangle = pixels
     where
 
         pixels :: [Pixel_]
@@ -146,7 +261,7 @@ getColorsInTriangle image (y', x') triangle = pixels
 
 
 getTriangleAverageRGB :: Image_ -> (Int, Int) -> S.Set (P2 Int)  -> C.Colour Double
-getTriangleAverageRGB image (y', x') triangle = blendEqually $ getColorsInTriangle image (y', x') triangle
+getTriangleAverageRGB image (x', y') triangle = blendEqually $ getColorsInTriangle image (x', y') triangle
 
 
 ptsBtween :: LineMXB -> LineMXB -> [(Int, Int)]