draw cusps in debug mode

src/app/MetaBrush/Render/Document.hs

@@ -72,7 +72,7 @@ import Math.Bezier.Spline
   , catMaybesSpline
   )
 import Math.Bezier.Stroke
-  ( CachedStroke(..), invalidateCache
+  ( Cusp(..), CachedStroke(..), invalidateCache
   , computeStrokeOutline
   )
 import Math.Interval
@@ -120,7 +120,7 @@ import MetaBrush.GTK.Util
 
 data Renders a
   = Renders
-  { renderStrokes, renderPath, renderBrushes
+  { renderStrokes, renderPath, renderDebug, renderBrushes
   , renderCLines, renderCPts, renderPPts :: a
   }
   deriving stock ( Show, Functor, Foldable, Traversable, Generic, Generic1 )
@@ -256,6 +256,7 @@ data StrokeRenderData where
                                     ( SplinePts Closed )
                                     ( SplinePts Closed, SplinePts Closed )
                                   , Seq FitPoint
+                                  , [ Cusp ]
                                   )
         , strokeBrushFunction :: pointParams -> SplinePts Closed
         }
@@ -299,14 +300,14 @@ strokeRenderData fitParams
                     let embedUsedParams = inject brush_defaults
 
                   -- Compute the outline using the brush function.
-                    ( outline, fitPts ) <-
+                    ( outline, fitPts, cusps ) <-
                       computeStrokeOutline @clo fitParams
                         ( toUsedParams . brushParams ) embedUsedParams brushFn
                         spline
                     pure $
                       StrokeWithOutlineRenderData
                         { strokeDataSpline    = spline
-                        , strokeOutlineData   = ( outline, fitPts )
+                        , strokeOutlineData   = ( outline, fitPts, cusps )
                         , strokeBrushFunction = fun @Double ( brushFn @Point proxy# id )
                                               . embedUsedParams
                                               . toUsedParams
@@ -324,13 +325,17 @@ renderStroke
 renderStroke cols@( Colours { brush } ) mbHoverContext mode rdrMode debug zoom = \case
   StrokeRenderData { strokeDataSpline } ->
     renderStrokeSpline cols mode rdrMode mbHoverContext zoom ( const ( pure () ) ) strokeDataSpline
-  StrokeWithOutlineRenderData strokeDataSpline strokeOutlineData strokeBrushFunction ->
+  StrokeWithOutlineRenderData strokeDataSpline ( strokeOutlineData, fitPts, cusps ) strokeBrushFunction ->
     renderStrokeSpline cols mode rdrMode mbHoverContext zoom
       ( when ( mode == BrushMode )
       . renderBrushShape ( cols { path = brush } ) mbHoverContext ( 2 * zoom ) strokeBrushFunction
       )
       strokeDataSpline
-    *> Compose blank { renderStrokes = drawOutline cols debug zoom strokeOutlineData }
+    *> Compose blank
+        { renderStrokes = drawOutline cols debug zoom strokeOutlineData
+        , renderDebug   =
+            when debug $ drawDebugInfo cols zoom ( fitPts, cusps )
+        }
 
 -- | Render a sequence of stroke points.
 --
@@ -550,13 +555,13 @@ drawCubicBezier ( Colours { path } ) zoom
 
 drawOutline
   :: Colours -> Bool -> Double
-  -> ( Either ( SplinePts Closed ) ( SplinePts Closed, SplinePts Closed ), Seq FitPoint )
+  -> Either ( SplinePts Closed ) ( SplinePts Closed, SplinePts Closed )
   -> Cairo.Render ()
-drawOutline cols@( Colours {..} ) debug zoom strokeData = do
+drawOutline ( Colours {..} ) debug zoom strokeData = do
   Cairo.save
   withRGBA brushStroke Cairo.setSourceRGBA
   case strokeData of
-    ( Left outline, fitPts ) -> do
+    Left outline -> do
       makeOutline outline
       case debug of
         False -> Cairo.fill
@@ -565,8 +570,7 @@ drawOutline cols@( Colours {..} ) debug zoom strokeData = do
           Cairo.setSourceRGBA 0 0 0 0.75
           Cairo.setLineWidth ( 2 / zoom )
           Cairo.stroke
-          ( `evalStateT` 0 ) $ traverse_ ( drawFitPoint cols zoom ) fitPts
+    Right ( fwd, bwd ) -> do
-    ( Right ( fwd, bwd ), fitPts ) -> do
       makeOutline fwd
       makeOutline bwd
       case debug of
@@ -576,7 +580,6 @@ drawOutline cols@( Colours {..} ) debug zoom strokeData = do
           Cairo.setSourceRGBA 0 0 0 0.75
           Cairo.setLineWidth ( 2 / zoom )
           Cairo.stroke
-          ( `evalStateT` 0 ) $ traverse_ ( drawFitPoint cols zoom ) fitPts
   Cairo.restore
   where
     makeOutline :: SplinePts Closed -> Cairo.Render ()
@@ -600,6 +603,14 @@ drawOutline cols@( Colours {..} ) debug zoom strokeData = do
         let ℝ2 x3 y3 = fromNextPoint start mp3
         in Cairo.curveTo x1 y1 x2 y2 x3 y3
 
+drawDebugInfo :: Colours -> Double
+              -> ( Seq FitPoint, [ Cusp ] )
+              -> Cairo.Render ()
+drawDebugInfo cols zoom ( fitPts, cusps ) = do
+  Cairo.setLineWidth ( 2 / zoom )
+  ( `evalStateT` 0 ) $ traverse_ ( drawFitPoint cols zoom ) fitPts
+  for_ cusps ( drawCusp cols zoom )
+
 drawFitPoint :: Colours -> Double -> FitPoint -> StateT Double Cairo.Render ()
 drawFitPoint _ zoom ( FitPoint { fitPoint = ℝ2 x y } ) = do
 
@@ -635,6 +646,27 @@ drawFitPoint _ zoom ( FitTangent { fitPoint = ℝ2 x y, fitTangent = V2 tx ty }
     Cairo.fill
     Cairo.restore
 
+drawCusp :: Colours -> Double -> Cusp -> Cairo.Render ()
+drawCusp _ zoom
+  ( Cusp { cuspPathCoords   = ℝ2 px py
+         , cuspStrokeCoords = ℝ2 cx cy } ) = do
+
+    -- Draw a white circle on the path at the cusp point.
+    Cairo.save
+    Cairo.translate px py
+    Cairo.arc 0 0 ( 4 / zoom ) 0 ( 2 * pi )
+    Cairo.setSourceRGBA 1.0 1.0 1.0 1.0
+    Cairo.stroke
+    Cairo.restore
+
+    -- Draw a black circle on the outline at the cusp point.
+    Cairo.save
+    Cairo.translate cx cy
+    Cairo.arc 0 0 ( 4 / zoom ) 0 ( 2 * pi )
+    Cairo.setSourceRGBA 0 0 0 1.0
+    Cairo.stroke
+    Cairo.restore
+
 drawSelectionRectangle :: Colours -> Double -> ℝ 2 -> ℝ 2 -> Cairo.Render ()
 drawSelectionRectangle ( Colours {..} ) zoom ( ℝ2 x0 y0 ) ( ℝ2 x1 y1 ) = do
 

src/splines/Math/Bezier/Stroke.hs

@@ -4,7 +4,7 @@
 {-# LANGUAGE UndecidableInstances  #-}
 
 module Math.Bezier.Stroke
-  ( Offset(..)
+  ( Offset(..), Cusp(..)
   , CachedStroke(..), discardCache, invalidateCache
   , computeStrokeOutline, joinWithBrush
   , withTangent
@@ -49,7 +49,7 @@ import GHC.Exts
 import GHC.STRef
   ( STRef(..), readSTRef, writeSTRef )
 import GHC.Generics
-  ( Generic, Generic1 )
+  ( Generic, Generic1, Generically(..) )
 import GHC.TypeNats
   ( type (-) )
 
@@ -158,13 +158,21 @@ data TwoSided a
   deriving stock    ( Show, Generic, Generic1, Functor, Foldable, Traversable )
   deriving anyclass ( NFData, NFData1 )
 
-type OutlineData = TwoSided ( SplinePts Open, Seq FitPoint )
+data OutlineData =
-instance Semigroup OutlineData where
+  OutlineData
+    { outline :: !( TwoSided ( SplinePts Open, Seq FitPoint ) )
+    , cusps   :: ![ Cusp ] }
+  deriving stock Generic
+  deriving ( Semigroup, Monoid )
+    via Generically OutlineData
+  deriving anyclass NFData
+
+instance Semigroup ( TwoSided ( SplinePts Open, Seq FitPoint ) ) where
   TwoSided ( fwdSpline1, fwdPts1 ) ( bwdSpline1, bwdPts1 ) <> TwoSided ( fwdSpline2, fwdPts2 ) ( bwdSpline2, bwdPts2 ) =
     TwoSided
       ( fwdSpline1 <> fwdSpline2, fwdPts1 <> fwdPts2 )
       ( bwdSpline2 <> bwdSpline1, bwdPts2 <> bwdPts1 )
-instance Monoid OutlineData where
+instance Monoid ( TwoSided ( SplinePts Open, Seq FitPoint ) ) where
   mempty = TwoSided empt empt
     where
       empt :: ( SplinePts Open, Seq FitPoint )
@@ -193,8 +201,26 @@ coords = view typed
 
 --------------------------------------------------------------------------------
 
+-- | Forward and backward outlines.
 type OutlineFn = ℝ 1 -> ( ( ℝ 2, T ( ℝ 2 ) ), ( ℝ 2, T ( ℝ 2 ) ) )
 
+data Cusp
+  = Cusp
+  { cuspParameters   :: !( ℝ 2 )
+    -- ^ @(t,s)@ parameter values of the cusp
+  , cuspPathCoords   :: !( ℝ 2 )
+    -- ^ path point coordinates
+  , cuspStrokeCoords :: !( ℝ 2 )
+    -- ^ brush stroke point coordinates
+  }
+  deriving stock Generic
+  deriving anyclass NFData
+
+data OutlineInfo =
+  OutlineInfo
+    { outlineFn :: OutlineFn
+    , outlineDefiniteCusps, outlinePotentialCusps :: [ Cusp ] }
+
 computeStrokeOutline ::
   forall ( clo :: SplineType ) usedParams brushParams crvData ptData s
   .  ( KnownSplineType clo
@@ -232,6 +258,7 @@ computeStrokeOutline ::
   -> ST s
         ( Either ( SplinePts Closed ) ( SplinePts Closed, SplinePts Closed )
         , Seq FitPoint
+        , [ Cusp ]
         )
 computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline { splineStart = spt0 } ) = case ssplineType @clo of
   -- Open brush path with at least one segment.
@@ -245,8 +272,8 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
         endPt :: ptData
         endPt = openCurveEnd lastCurve
         startTgtFwd, startTgtBwd, endTgtFwd, endTgtBwd :: T ( ℝ 2 )
-        ( ( _, startTgtFwd), ( _, startTgtBwd ) ) = firstOutlineFn $ ℝ1 0
+        ( ( _, startTgtFwd), ( _, startTgtBwd ) ) = outlineFn firstOutlineFn $ ℝ1 0
-        ( ( _, endTgtFwd  ), ( _, endTgtBwd   ) ) =  lastOutlineFn $ ℝ1 1
+        ( ( _, endTgtFwd  ), ( _, endTgtBwd   ) ) = outlineFn  lastOutlineFn $ ℝ1 1
         startBrush, endBrush :: SplinePts Closed
         startBrush = brushShape spt0
         endBrush   = brushShape endPt
@@ -286,10 +313,14 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
           $ joinWithBrush endBrush endTgtBwd endTgtFwd
 
     -> do
-      TwoSided ( fwdPts, fwdFits ) ( bwdPts, bwdFits ) <- updateSpline ( startTgt, startTgtFwd, startTgtBwd )
+      OutlineData
+        ( TwoSided ( fwdPts, fwdFits ) ( bwdPts, bwdFits ) )
+        cusps
+          <- updateSpline ( startTgt, startTgtFwd, startTgtBwd )
       pure
         ( Left ( adjustSplineType @Closed $ startCap <> fwdPts <> endCap <> bwdPts )
         , fwdFits <> bwdFits
+        , cusps
         )
   -- Closed brush path with at least one segment.
   -- Add forward and backward caps at the start.
@@ -305,33 +336,39 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
         endTgt = case prevCurves of
           Empty          -> endTangent spt0 spt0                      lastCurve
           _ :|> lastPrev -> endTangent spt0 ( openCurveEnd lastPrev ) lastCurve
-        ( ( _, startTgtFwd), ( _, startTgtBwd ) ) = firstOutlineFn $ ℝ1 0
+        ( ( _, startTgtFwd), ( _, startTgtBwd ) ) = outlineFn firstOutlineFn $ ℝ1 0
-        ( ( _, endTgtFwd  ), ( _, endTgtBwd   ) ) =  lastOutlineFn $ ℝ1 1
+        ( ( _, endTgtFwd  ), ( _, endTgtBwd   ) ) = outlineFn  lastOutlineFn $ ℝ1 1
         fwdStartCap, bwdStartCap :: SplinePts Open
-        TwoSided fwdStartCap bwdStartCap
+        cusps :: [ Cusp ]
-          = fmap fst . snd . runWriter
+        OutlineData ( fmap fst -> TwoSided fwdStartCap bwdStartCap ) cusps
+          = snd . runWriter
           $ tellBrushJoin ( endTgt, endTgtFwd, endTgtBwd ) spt0 ( startTgt, startTgtFwd, startTgtBwd )
     -> do
-      TwoSided ( fwdPts, fwdFits ) ( bwdPts, bwdFits ) <- updateSpline ( endTgt, endTgtFwd, endTgtBwd )
+      OutlineData
+        ( TwoSided ( fwdPts, fwdFits ) ( bwdPts, bwdFits ) )
+        cusps
+          <- updateSpline ( endTgt, endTgtFwd, endTgtBwd )
       pure
         ( Right ( adjustSplineType @Closed ( fwdStartCap <> fwdPts ), adjustSplineType @Closed ( bwdPts <> bwdStartCap ) )
         , fwdFits <> bwdFits
+        , cusps
         )
   -- Single point.
   _ ->
     pure
       ( Left $ fmap ( T ( coords spt0 ) • ) ( brushShape spt0 )
       , Empty
+      , []
       )
   where
 
-    outlineFns :: Seq OutlineFn
+    outlineFns :: Seq OutlineInfo
     outlineFns = go spt0 ( openCurves $ splineCurves ( adjustSplineType @Open spline ) )
       where
         go
           :: ptData
           -> Seq ( Curve Open crvData ptData )
-          -> Seq OutlineFn
+          -> Seq OutlineInfo
         go _ Empty = Empty
         go p0 ( crv :<| crvs ) =
           outlineFunction ptParams toBrushParams brushFn p0 crv :<| go ( openCurveEnd crv ) crvs
@@ -347,13 +384,13 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
           ( \ ptData curve -> do
             ( prevTgt, prev_tgtFwd, prev_tgtBwd ) <- get
             let
-              fwdBwd :: OutlineFn
+              fwdBwd :: OutlineInfo
               fwdBwd = outlineFunction ptParams toBrushParams brushFn ptData curve
               tgt, next_tgt, tgtFwd, next_tgtFwd, tgtBwd, next_tgtBwd :: T ( ℝ 2 )
               tgt      = startTangent spt0 ptData curve
               next_tgt =   endTangent spt0 ptData curve
-              ( ( _, tgtFwd      ), ( _, tgtBwd      ) ) = fwdBwd $ ℝ1 0
+              ( ( _, tgtFwd      ), ( _, tgtBwd      ) ) = outlineFn fwdBwd $ ℝ1 0
-              ( ( _, next_tgtFwd ), ( _, next_tgtBwd ) ) = fwdBwd $ ℝ1 1
+              ( ( _, next_tgtFwd ), ( _, next_tgtBwd ) ) = outlineFn fwdBwd $ ℝ1 1
             lift $ tellBrushJoin ( prevTgt, prev_tgtFwd, tgtBwd ) ptData ( tgt, tgtFwd, prev_tgtBwd )
             lift $ updateCurveData ( curveData curve ) fwdBwd
             put ( next_tgt, next_tgtFwd, next_tgtBwd )
@@ -363,7 +400,7 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
 
     updateCurveData
       :: crvData
-      -> OutlineFn
+      -> OutlineInfo
       -> WriterT OutlineData ( ST s ) ()
     updateCurveData ( view ( typed @( CachedStroke s ) ) -> CachedStroke { cachedStrokeRef } ) fwdBwd = do
       mbOutline <- lift ( readSTRef cachedStrokeRef )
@@ -375,13 +412,21 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
         -- No cached fit: compute the fit anew.
         Nothing -> do
           let
+            cusps :: [ Cusp ]
+            cusps = outlineDefiniteCusps  fwdBwd
+                 ++ outlinePotentialCusps fwdBwd
+            -- SLD TODO: here!!!
             fwdData, bwdData :: ( SplinePts Open, Seq FitPoint )
             ( fwdData, bwdData ) =
-              ( fitSpline fitParams ( fst . fwdBwd ), fitSpline fitParams ( snd . fwdBwd ) )
+              ( fitSpline fitParams ( fst . outlineFn fwdBwd )
+              , fitSpline fitParams ( snd . outlineFn fwdBwd ) )
                 `Strats.using`
                   ( Strats.parTuple2 Strats.rdeepseq Strats.rdeepseq )
             outlineData :: OutlineData
-            outlineData = TwoSided fwdData ( bimap reverseSpline Seq.reverse bwdData )
+            outlineData =
+              OutlineData
+                ( TwoSided fwdData ( bimap reverseSpline Seq.reverse bwdData ) )
+                cusps
           outlineData `deepseq` tell outlineData
           lift $ writeSTRef cachedStrokeRef ( Just outlineData )
 
@@ -394,7 +439,9 @@ computeStrokeOutline fitParams ptParams toBrushParams brushFn spline@( Spline {
       -> ( T ( ℝ 2 ), T ( ℝ 2 ), T ( ℝ 2 ) )
       -> WriterT OutlineData m ()
     tellBrushJoin ( prevTgt, prevTgtFwd, prevTgtBwd ) sp0 ( tgt, tgtFwd, tgtBwd ) =
-      tell $ TwoSided ( fwdJoin, Empty ) ( bwdJoin, Empty )
+      tell $ OutlineData
+               ( TwoSided ( fwdJoin, Empty ) ( bwdJoin, Empty ) )
+               []
       where
         ptOffset :: T ( ℝ 2 )
         ptOffset = ℝ2 0 0 --> coords sp0
@@ -477,7 +524,7 @@ outlineFunction
       )
   -> ptData
   -> Curve Open crvData ptData
-  -> OutlineFn
+  -> OutlineInfo
 outlineFunction ptParams toBrushParams brushFromParams sp0 crv =
   let
     pathAndUsedParams :: forall i k arr
@@ -511,6 +558,16 @@ outlineFunction ptParams toBrushParams brushFromParams sp0 crv =
     path :: C ( ExtentOrder 'Point ) ( ℝ 1 ) ( ℝ 2 )
     ( path, usedParams ) = pathAndUsedParams @Point id
 
+    curves :: ℝ 1 -> Seq ( ℝ 1 -> StrokeDatum Point )
+    curves =
+      brushStrokeData @Point @( ExtentOrder 'Point ) @brushParams
+        path
+        ( chainRule @Double @( ExtentOrder 'Point )
+            usedParams
+            ( linear toBrushParams )
+        )
+        ( brushFromParams @Point proxy# id )
+
     curvesI :: 𝕀ℝ 1 -> Seq ( 𝕀ℝ 1 -> StrokeDatum 'Interval )
     curvesI = brushStrokeData @'Interval @( ExtentOrder 'Interval ) @brushParams
                pathI
@@ -526,21 +583,10 @@ outlineFunction ptParams toBrushParams brushFromParams sp0 crv =
 
     fwdBwd :: OutlineFn
     fwdBwd t
-      = solveEnvelopeEquations t path_t path'_t ( fwdOffset, bwdOffset ) curves
+      = solveEnvelopeEquations t path_t path'_t ( fwdOffset, bwdOffset )
-    --  = ( ( offset fwdOffset • path_t,       path'_t )
+          ( curves t )
-    --    , ( offset bwdOffset • path_t, -1 *^ path'_t ) )
       where
 
-        curves :: Seq ( ℝ 1 -> StrokeDatum Point )
-        curves = brushStrokeData @Point @( ExtentOrder 'Point ) @brushParams
-                   path
-                   ( chainRule @Double @( ExtentOrder 'Point )
-                       usedParams
-                       ( linear toBrushParams )
-                   )
-                   ( brushFromParams @Point proxy# id )
-                   t
-
         fwdOffset = withTangent         path'_t   brush_t
         bwdOffset = withTangent ( -1 *^ path'_t ) brush_t
 
@@ -563,7 +609,10 @@ outlineFunction ptParams toBrushParams brushFromParams sp0 crv =
      --     , "unknown:"
      --     , if null newtDunno then "[]" else unlines $ map show newtDunno ]
      --   ) $
-     fwdBwd
+     OutlineInfo
+       { outlineFn = fwdBwd
+       , outlineDefiniteCusps  = map ( cuspCoords curves ) newtSols
+       , outlinePotentialCusps = map ( cuspCoords curves ) newtDunno }
 
 -----------------------------------
 -- Various utility functions
@@ -1117,6 +1166,26 @@ extendedRecip x@( 𝕀 lo hi )
   | otherwise
   = [ recip ( 𝕀 lo 0 ), recip ( 𝕀 0 hi ) ]
 
+-- | Computes the brush stroke coordinates of a cusp from
+-- the @(t,s)@ parameter values.
+cuspCoords :: ( ℝ 1 -> Seq ( ℝ 1 -> StrokeDatum 'Point ) )
+           -> ( 𝕀ℝ 1, Int, 𝕀ℝ 1 )
+           -> Cusp
+cuspCoords eqs ( 𝕀 ( ℝ1 t_lo ) ( ℝ1 t_hi ), i, 𝕀 ( ℝ1 s_lo ) ( ℝ1 s_hi ) )
+  | StrokeDatum
+    { dpath   = D21 { _D21_v = path  }
+    , dstroke = D22 { _D22_v = stroke } }
+      <- ( eqs ( ℝ1 t_mid ) `Seq.index` i ) ( ℝ1 s_mid )
+  = Cusp
+  { cuspParameters   = ℝ2 t_mid s_mid
+  , cuspPathCoords   = path
+  , cuspStrokeCoords = stroke
+  }
+  where
+    t_mid = 0.5 * ( t_lo + t_hi )
+    s_mid = 0.5 * ( t_lo + t_hi )
+
+
 -- | Interval Newton method with Gauss–Seidel step for inversion
 -- of the interval Jacobian.
 --

src/splines/Math/Interval/Internal.hs

@@ -16,6 +16,10 @@ import Data.Monoid
 import GHC.Show
   ( showCommaSpace )
 
+-- deepseq
+import Control.DeepSeq
+  ( NFData )
+
 -- rounded-hw
 import Numeric.Rounded.Hardware
   ( Rounded(..) )
@@ -40,6 +44,7 @@ import Math.Ring
 
 newtype 𝕀 a = MkI { ival :: Interval a }
   deriving newtype ( Prelude.Num, Prelude.Fractional, Prelude.Floating )
+  deriving newtype NFData
 
 instance Eq a => Eq ( 𝕀 a ) where
   𝕀 a b == 𝕀 c d =