typography-geometry 1.0 → 1.0.0
raw patch · 9 files changed
+1303/−1302 lines, 9 files
Files
- Graphics/Typography.lhs +0/−89
- Graphics/Typography/Approximation.lhs +0/−170
- Graphics/Typography/Bezier.lhs +0/−749
- Graphics/Typography/Geometry.lhs +89/−0
- Graphics/Typography/Geometry/Approximation.lhs +170/−0
- Graphics/Typography/Geometry/Bezier.lhs +749/−0
- Graphics/Typography/Geometry/Outlines.lhs +291/−0
- Graphics/Typography/Outlines.lhs +0/−291
- typography-geometry.cabal +4/−3
− Graphics/Typography.lhs
@@ -1,89 +0,0 @@-\begin{code}-{-# OPTIONS -XFlexibleInstances -XNamedFieldPuns #-}--- | This module contains basic tools for geometric types and functions.-module Graphics.Typography (Matrix2(..),- inverse,rotation,- Geometric(..),- leftMost,rightMost,topMost,bottomMost)- where--import Algebra.Polynomials.Numerical---- | The type of the transformation matrices used in all geometrical applications.-data Matrix2 a=- -- | The application of @Matrix2 a b c d@ to vector @(x,y)@ should be- -- @(ax+by,cx+dy)@.- Matrix2 a a a a deriving (Show, Read, Eq)---- | Inverses an inversible matrix. If it is not inversible,--- The behaviour is undefined.-inverse::(Fractional a, Num a)=>Matrix2 a->Matrix2 a-inverse (Matrix2 a b c d)=- let det=a*d-c*b in- Matrix2 (d/det) (-b/det) (-c/det) (a/det)-- --instance Num a=>Num (Matrix2 a) where- (+) (Matrix2 a b c d) (Matrix2 e f g h)=- Matrix2 (a+e) (b+f) (c+g) (d+h)- (*) (Matrix2 a b c d) (Matrix2 e f g h)=- Matrix2 (a*e+b*g) (a*f+b*h) (c*e+d*g) (c*f+d*h)- fromInteger a=Matrix2 (fromInteger a) 0 0 (fromInteger a)- abs=undefined- signum=undefined--instance Intervalize Matrix2 where- intervalize (Matrix2 a b c d)=- Matrix2 (interval a) (interval b) (interval c) (interval d)-- intersects (Matrix2 a b c d) (Matrix2 a' b' c' d')=- (intersectsd a a') &&- (intersectsd b b') &&- (intersectsd c c') &&- (intersectsd d d')- --- | A class for applying geometric applications to objects-class Geometric g where- translate::Double->Double->g->g- apply::Matrix2 Double->g->g---- | The matrix of a rotation-rotation::Floating a=>a->Matrix2 a-rotation theta=- let ct=cos theta- st=sin theta- in- Matrix2 ct (-st) st ct--instance Geometric g=>Geometric [g] where- - translate x y cur=map (translate x y) cur- apply m cur=map (apply m) cur- ---- | @'leftMost' a b@ is the leftmost point between @a@ and @b@.-leftMost::(Double,Double)->(Double,Double)->(Double,Double)-leftMost u@(a,_) v@(b,_)- | a<b = u- | otherwise = v--- | @'rightMost' a b@ is the rightmost point between @a@ and @b@.-rightMost::(Double,Double)->(Double,Double)->(Double,Double)-rightMost u@(a,_) v@(b,_)- | a<b = v- | otherwise = u--- | @'bottomMost' a b@ is the lower point between @a@ and @b@.-bottomMost::(Double,Double)->(Double,Double)->(Double,Double)-bottomMost u@(_,a) v@(_,b)- | a<b = u- | otherwise = v--- | @'topMost' a b@ is the upper point between @a@ and @b@.-topMost::(Double,Double)->(Double,Double)->(Double,Double)-topMost u@(_,a) v@(_,b)- | a<b = v- | otherwise = u-----\end{code}
− Graphics/Typography/Approximation.lhs
@@ -1,170 +0,0 @@-\begin{code}-{-# OPTIONS -XRecordWildCards -XNamedFieldPuns #-}--- | This module contains the function to approximate a list of curves with--- degree 3 Bezier curves, using a least squares method.--module Graphics.Typography.Approximation(approximate) where--import qualified Data.Vector.Unboxed as UV-import Graphics.Typography.Bezier-import Graphics.Typography-import Algebra.Polynomials.Bernstein--import Algebra.Polynomials.Numerical--- import Debug.Trace-rnd::Interval->Double-rnd (Interval a b)=(a+b)/2----- | Approximates a list of 'Curves' with a list of degree 3 Bernstein curves.-approximate::[Curve]->[(Bernsteinp Int Double,Bernsteinp Int Double)]-approximate []=[]-approximate l0@(h0:_)= -- traceShow "starting" $- let approx::Double->Double->[Curve]->[(Bernsteinp Int Double,Bernsteinp Int Double)]- approx _ _ []=[]- approx x0 y0 (cc@(Circle {..}):s)= -- traceShow "circle" $- let theta=abs $ t1-t0 in- if theta <= pi/2 then- let x0_=cos $ theta/2- y0_=sin $ theta/2- x1_=(4-x0_)/3- y1_=(1-x0_)*(3-x0_)/(3*y0_)- - c0=cos $! theta/2+t0- s0=sin $! theta/2+t0- - - px0=c0*x0_ - s0*y0_- py0=s0*x0_ + c0*y0_- - px1=c0*x1_ - s0*y1_- py1=s0*x1_ + c0*y1_- - px2=c0*x1_ + s0*y1_- py2=s0*x1_ - c0*y1_- - -- px3=c0*x0_ + s0*y0_- -- py3=s0*x0_ - c0*y0_- - x1=cx0+(a*px0+b*py0)- y1=cy0+(c*px0+d*py0)- - (Matrix2 a b c d)=matrix- x=Bernsteinp 4 $ UV.fromList- [ x0, -- cx0+(a*px3+b*py3),- cx0+(a*px2+b*py2),- cx0+(a*px1+b*py1),- x1]- - y=Bernsteinp 4 $ UV.fromList- [ y0, -- cy0+(c*px3+d*py3),- cy0+(c*px2+d*py2),- cy0+(c*px1+d*py1),- y1 ]- in- (x,y):(approx x1 y1 s)- - else- let t1'=(t1+t0)/2 in- approx x0 y0 $ (cc { t1=t1' }):(cc { t0=t1' }):s-{-- approx x0 y0 (h@(Bezier{}):s)=- -- incorrect !- (restriction (cx h) (t0 h) (t1 h),- restriction (cy h) (t0 h) (t1 h)):- (approx (UV.last $ coefs $ cx h)- (UV.last $ coefs $ cy h) s)--}- -- Ce qui suit est une methode de moindres carres- approx x0 y0 (off_:s)= -- traceShow ("offset") $- -- On commence par chercher les points ou la derivee de la norme- -- de la tangente est maximale. C'est la qu'on va couper s'il y- -- a un probleme.- let bx=restriction (cx off_) (t0 off_) (t1 off_)- by=restriction (cy off_) (t0 off_) (t1 off_)- off=off_ { cx=bx,cy=by,t0=0,t1=1 }- ibx=elevate (bounds by-bounds bx) $ intervalize bx- iby=elevate (bounds bx-bounds by) $ intervalize by- points=- let np=10 in- map (\x->(x/np,x/np)) [0..np]- -- Ensuite, moindres carres standard, comme dans Hoschek 88.- - vx0=ibx?1-ibx?0- vy0=iby?1-iby?0- vx1=ibx?(bounds ibx-2)-ibx?(bounds ibx-1)- vy1=iby?(bounds iby-2)-iby?(bounds iby-1)- - (wx0,wy0)=evalCurve off 0- (wx1,wy1)=evalCurve off 1-- wx=Bernsteinp 4 $ UV.fromList [wx0,wx0,wx1,wx1] :: Bernsteinp Int Interval- wy=Bernsteinp 4 $ UV.fromList [wy0,wy0,wy1,wy1] :: Bernsteinp Int Interval-- bern1=Bernsteinp 4 $ UV.fromList [0,1,0,0] :: Bernsteinp Int Interval- bern2=Bernsteinp 4 $ UV.fromList [0,0,1,0] :: Bernsteinp Int Interval-- sumAll a b c d x1 y1 ((h1,h2):ss)=- - let h=Interval h1 h2- (xi,yi)=evalCurve off h- - b1=eval bern1 h- b2=eval bern2 h- - a'=a + (vx0*vx0+vy0*vy0)*b1*b1- b'=b + (vx0*vx1 + vy0*vy1)*b1*b2- c'=c + (vx0*vx1 + vy0*vy1)*b1*b2- d'=d + (vx1*vx1+vy1*vy1)*b2*b2- - dx=xi-(eval wx h)- dy=yi-(eval wy h)- - x1'=x1 + (vx0*dx + vy0*dy)*b1- y1'=y1 + (vx1*dx + vy1*dy)*b2- in- sumAll a' b' c' d' x1' y1' ss- - sumAll a b c d x1 y1 []=(a,b,c,d,x1,y1)- - (ra,rb,rc,rd,rx1,ry1)=sumAll 0 0 0 0 0 0 points- - (Matrix2 a_ b_ c_ d_)=inverse $ Matrix2 ra rb rc rd- lambda1=a_*rx1+b_*ry1- lambda2=c_*rx1+d_*ry1- - -- On a la courbe optimale. Il faut chercher ou on va couper, maintenant- xap=Bernsteinp 4 $ UV.fromList [wx0,- wx0+lambda1*vx0,- wx1+lambda2*vx1,- wx1]- yap=Bernsteinp 4 $ UV.fromList [wy0,- wy0+lambda1*vy0,- wy1+lambda2*vy1,- wy1]- (err,arg)=foldl (\m (h1,h2)->- let (xi,yi)=evalCurve off (Interval h1 h2)- xj=eval xap (Interval h1 h2)- yj=eval yap (Interval h1 h2)- in- max m (iup $ abs $ (xi-xj)*(xi-xj)+(yi-yj)*(yi-yj), (h1+h2)/2))- (0,0) points- in- if err<=0.01 then- (desintervalize xap,desintervalize yap):(approx (rnd wx1) (rnd wy1) s)- else- approx x0 y0 $- (off { cx=restriction (cx off) 0 arg,- cy=restriction (cy off) 0 arg }):- (off { cx=restriction (cx off) arg 1,- cy=restriction (cy off) arg 1 }):s- - (x0h,y0h)=evalCurve h0 $ Interval (t0 h0) (t0 h0)- - in- approx (rnd x0h) (rnd y0h) l0- -desintervalize::(Bernsteinp a Interval)->(Bernsteinp a Double)-desintervalize b=b { coefs=UV.map rnd $ coefs b}- -\end{code}
− Graphics/Typography/Bezier.lhs
@@ -1,749 +0,0 @@-\documentclass{article}-%include lhs2TeX.fmt-\begin{document}-\begin{code}-{-# OPTIONS -XUnboxedTuples -XBangPatterns -XNamedFieldPuns -XRecordWildCards -XMagicHash -cpp #-}--- | This module contains the basic functions for manipulating Bezier curves. It is heavily--- based on the book by N. M. Patrikalakis and T. Maekawa, Shape Interrogation for Computer--- Aided Design and Manufacturing.--module Graphics.Typography.Bezier (- Curve(..),line,bezier3,- offset,- inter,- evalCurve,distance,- left,bottom,right,top) where--import Algebra.Polynomials.Bernstein-import Algebra.Polynomials.Numerical-import Graphics.Typography-import qualified Data.Vector as V-import qualified Data.Vector.Unboxed as UV-import Data.List (partition,sort)---- | The type for representing all types of curves.-data Curve=- Bezier { cx::Bernsteinp Int Double,- cy::Bernsteinp Int Double, - t0::Double,- t1::Double }- - | Offset { cx::Bernsteinp Int Double,- cy::Bernsteinp Int Double,- t0::Double,- t1::Double,- matrix::Matrix2 Double- }- | Circle { cx0::Double,- cy0::Double,- t0::Double,- t1::Double,- matrix::Matrix2 Double- }- deriving (Show)----- | The basic constructor for lines : a line is a degree 1 Bezier curve-line::Double->Double->Double->Double->Curve-line px py px' py'=Bezier { cx=Bernsteinp 2 $ UV.fromList [px,px'], - cy=Bernsteinp 2 $ UV.fromList [py,py'],- t0=0,t1=1 }---- | A shortcut to define degree 3 Bezier curves from points. If the control--- points are @a,b,c,d@, the function should be called with--- @'bezier3' xa ya xb yb xc yc xd yd@.-bezier3::Double->Double->Double->Double->Double->Double->Double->Double->Curve-bezier3 px0 py0 px1 py1 px2 py2 px3 py3=- Bezier { cx=Bernsteinp 4 $ UV.fromList [px0,px1,px2,px3],- cy=Bernsteinp 4 $ UV.fromList [py0,py1,py2,py3],- t0=0,t1=1 }--instance Geometric Curve where- translate x y cur@(Circle{cx0,cy0})= - cur { cx0=cx0+x,cy0=cy0+y }- translate x y cur=- cur { cx=(cx cur) { coefs=UV.map (+x) $ coefs $ cx cur},- cy=(cy cur) { coefs=UV.map (+y) $ coefs $ cy cur} }-- apply m0@(Matrix2 a b c d) cir@(Circle{cx0,cy0,matrix})=- cir { cx0=a*cx0+b*cy0, cy0=c*cx0+d*cy0, matrix=m0*matrix }- apply (Matrix2 a b c d) cur=- cur { cx=(scale a $ cx cur)+(scale b $ cy cur),- cy=(scale c $ cx cur)+(scale d $ cy cur) }------ | Gives the point corresponding to the given value of the parameter-evalCurve::Curve->Interval->(Interval,Interval)-evalCurve (Offset{..}) t=- let ix=intervalize cx- iy=intervalize cy- xt=eval ix t- yt=eval iy t- m@(Matrix2 a b c d)=intervalize matrix- (Matrix2 a_ b_ c_ d_)=inverse m- xt0'=eval (derivate ix) t- yt0'=eval (derivate iy) t- xt'=a_*xt0' + b_*yt0'- yt'=c_*xt0' + d_*yt0'- dd=sqrt $ xt'*xt' + yt'*yt'- in- (xt+(a*yt'-b*xt')/dd, yt+(c*yt'-d*xt')/dd)--evalCurve (Circle{..}) alpha=- let xx=cos alpha- yy=sin alpha- (Matrix2 a b c d)=intervalize matrix- in- (interval cx0+a*xx+b*yy, interval cy0+c*xx+d*yy)- -evalCurve (Bezier{..}) t=- let ix=intervalize cx- iy=intervalize cy- xx=eval ix t- yy=eval iy t- in- (xx,yy)- -data Topo=Dehors | SurLaLigne | Dedans deriving Eq---- | @'inter' c0 c1@ is a list of all possible points of intersection--- between curves @c0@ and @c1@ : if @(u,v,w,x)@ is returned by 'inter',--- then curve @c0@ may intersect with @c1@ between parameter values @u@--- and @v@, which corresponds to parameter values between @w@ and @x@ for--- @c1@. The implementation guarantees that all actual solutions are found,--- but possibly false solutions may also be returned.--inter::Curve->Curve->[((Double,Double,Double,Double))]-inter op@(Offset { cx=bxp_,cy=byp_,matrix=mp,t0=t0a,t1=t1a })- (Offset { cx=bxq_,cy=byq_,matrix=mq,t0=t0b,t1=t1b })=- - -- Attention : verifier si c'est la meme generatrice- let thrx=1e-5- solutions _ []=[]- solutions thr boxes@(_:_)=- let sol0=concatMap (solve thr (V.fromList [eq0,eq1,eq2,eq3])) boxes - - (correct,toRefine)=partition (\(u,v,_,_,_,_,_,_)->- let (xu,yu)=evalCurve op (Interval u u)- (xv,yv)=evalCurve op (Interval v v)- in- (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0- in- correct++(solutions (thr/2) toRefine)- in- map (\(u,v,w,x,_,_,_,_)->(u,v,w,x)) $ solutions 1e-2 $- [(t0a,t1a,t0b,t1b,0,1,0,1)::- (Double,Double,Double,Double,Double,Double,Double,Double)]-- where- - - imp@(Matrix2 ap bp cp dp)=intervalize mp- imq@(Matrix2 aq bq cq dq)=intervalize mq- (Matrix2 ap_ bp_ cp_ dp_)=inverse imp- (Matrix2 aq_ bq_ cq_ dq_)=inverse imq-- bxp=intervalize bxp_- byp=intervalize byp_- bxq=intervalize bxq_- byq=intervalize byq_- - bxp4=promote 1 bxp- byp4=promote 1 byp- bxq4=promote 2 bxq- byq4=promote 2 byq- -- bxp'=derivate bxp- byp'=derivate byp- bxq'=derivate bxq- byq'=derivate byq-- bXp'=promote 1 $ (scale ap_ bxp')+(scale bp_ byp')- bYp'=promote 1 $ (scale cp_ bxp')+(scale dp_ byp')-- bXq'=promote 2 $ (scale aq_ bxq')+(scale bq_ byq')- bYq'=promote 2 $ (scale cq_ bxq')+(scale dq_ byq')-- bomp@(Bernsteinp _ omegap)=(bXp'*bXp')+(bYp'*bYp') :: Bernsteinp (Int,Int,Int,Int) Interval- bomq@(Bernsteinp _ omegaq)=(bXq'*bXq')+(bYq'*bYq') :: Bernsteinp (Int,Int,Int,Int) Interval- - au=- let au_=sqrt $ UV.minimum $ UV.map ilow omegap in- max 0 $ fpred au_- bu=- let bu_=sqrt $ UV.maximum $ UV.map iup omegap in- fsucc bu_- av=- let av_=sqrt $ UV.minimum $ UV.map ilow omegaq in- max 0 $ fpred av_- bv=- let bv_=sqrt $ UV.maximum $ UV.map iup omegaq in- fsucc bv_- - alphau=Bernsteinp (1,1,2,1) $ UV.fromList [Interval au au,Interval bu bu]- alphav=Bernsteinp (1,1,1,2) $ UV.fromList [Interval av av,Interval bv bv]-- eq0=- ((bxp4*alphau*alphav) + (scale ap $ bYp'*alphav) - (scale bp $ bXp'*alphav)- -(bxq4*alphau*alphav) - (scale aq $ bYq'*alphau) + (scale bq $ bXq'*alphau))- eq1=- ((byp4*alphau*alphav) + (scale cp $ bYp'*alphav) - (scale dp $ bXp'*alphav)- -(byq4*alphau*alphav) - (scale cq $ bYq'*alphau) + (scale dq $ bXq'*alphau)) - eq2=bomp-(alphau*alphau)- eq3=bomq-(alphav*alphav)- - ---inter b@(Circle{}) a@(Offset{})=- map (\(i,j,k,l)->(k,l,i,j)) $ inter a b--inter o@(Offset { cx=bxp, cy=byp, matrix=mp })- cir@(Circle{cx0,cy0,matrix=mq})=-- let ix=intervalize bxp- iy=intervalize byp- m@(Matrix2 a b c d)=intervalize mp- (Matrix2 a_ b_ c_ d_)=inverse m- x'=derivate ix- y'=derivate iy- xx'=(scale a_ x')+(scale b_ y')- yy'=(scale c_ x')+(scale d_ y')- omega@(Bernsteinp _ omegap)=xx'*xx'+yy'*yy'- au=- let au_=sqrt $ UV.minimum $ UV.map ilow omegap in- max 0 $ fpred au_- bu=- let bu_=sqrt $ UV.maximum $ UV.map iup omegap in- fsucc bu_- - alphau=Bernsteinp (1,2) $ UV.fromList [Interval au au,Interval bu bu]- - lambda=(promote 1 omega) - alphau*alphau- -- Avant multiplication par M_C^-1- xx0=(promote 1 $ ix-(intervalize $ constant cx0))*alphau- +(promote 1 $ scale a yy'-scale b xx')- yy0=(promote 1 $ iy-(intervalize $ constant cy0))*alphau- +(promote 1 $ scale c yy'-scale d xx') :: Bernsteinp (Int,Int) Interval- - (Matrix2 ac_ bc_ cc_ dc_)=inverse $ intervalize mq- xx1=(scale ac_ xx0)+(scale bc_ yy0)- yy1=(scale cc_ xx0)+(scale dc_ yy0)- - eqc=xx1*xx1+yy1*yy1-alphau*alphau- - thrx=1e-5- - solutions _ []=[]- solutions thr boxes@(_:_)=- let sol0=concatMap (solve thr (V.fromList [eqc,lambda])) boxes - - (correct,toRefine)=partition (\(u,v,_,_)->- let (xu,yu)=evalCurve o (Interval u u)- (xv,yv)=evalCurve o (Interval v v)- in- (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0- in- correct++(solutions (thr/2) toRefine)- -- -- Removing false positives by computing the distance to the center of- -- the circle (this is quite fast).- - removeFalse cl0 (h@(_,v,_,_):h'@(u',_,_,_):s)=- let u''=(v+u')/2- (xu,yu)=evalCurve o (Interval u'' u'')- Interval dl du=distance xu yu cir- cl1- | du<1 = Dedans- | dl>1 = Dehors- | otherwise = SurLaLigne- in- if cl0/=cl1 then h:(removeFalse cl1 (h':s)) else- removeFalse cl1 (h':s)- removeFalse _ l=l- - initCl=- let (x0,y0)=evalCurve o (Interval (t0 o) (t0 o)) - Interval dl du=distance x0 y0 cir- in- if dl>1 then Dehors else if du<1 then Dedans else SurLaLigne- in- foldl (\l (u,v,_,_)->- let (Interval xl xu,Interval yl yu)=evalCurve o (Interval u v) in- case angle (Interval xl xu) (Interval yl yu) cir of- Just (Interval a0 a1)->- (u,v,a0,a1):l- Nothing->l- ) [] $ removeFalse initCl $ sort $ solutions 1e-2 [(t0 o,t1 o,0::Double,1::Double)]- - -inter a@(Circle{cx0=x0a,cy0=y0a,matrix=ma})- b@(Circle{cx0=x0b,cy0=y0b,matrix=mb})=- - if (intervalize ma)`intersects`(intervalize mb) && x0a==x0b && y0a==y0b then- let up ix@(Interval _ x_) tt0 tt1- | x_<tt0 =- up (ix+(2*interval pi)) tt0 tt1- | otherwise = down ix tt0 tt1- down ix@(Interval x_ x__) tt0 tt1- | x_>tt1 =- down (ix-(2*interval pi)) tt0 tt1- | x__<tt0 =- Nothing- | otherwise =- Just ix- - alpha=up (interval $ t0 a) (t0 b) (t1 b)- beta=up (interval $ t0 b) (t0 b) (t1 b)- in- - case (alpha,beta) of- (Just aa,Just ab)->- case (up aa (t0 a) (t1 a),- up ab (t0 a) (t1 a)) of- - (Just ba,Just bb)- | ilow aa<=iup ab -> [(ilow ba, iup bb,- ilow aa, iup ab)]- | otherwise->- case (up (interval $ t0 b) (t0 a) (t1 a),- up (interval $ t1 b) (t0 a) (t1 a)) of- (Just b0,Just b1)->- [(ilow b0,iup bb,- t0 b, iup ab),- (ilow ba,iup b1,- ilow aa, t1 b)]- _->[]- _->[]- _->[]-- else- let thr=1e-5- solutions=solve thr (V.fromList [eq0,eq1]) (fpred u0,fsucc v0,- fpred w0,fsucc x0)- in- foldl (\l (u,v,w,x)->- let alpha=angle (Interval u v) (Interval w x) a- beta=angle (Interval u v) (Interval w x) b- in- case alpha of- Just (Interval a0l a0u)->- case beta of- Just (Interval b0l b0u)->(a0l,a0u,b0l,b0u):l- _->l- _->l- ) [] solutions - where- - ima@(Matrix2 am bm cm dm)=intervalize ma- - maxa=max (iup $ abs am+abs bm) (iup $ abs cm+abs dm)- (u0,v0,w0,x0)=(x0a-maxa,x0a+maxa,y0a-maxa,y0a+maxa)- - -- x-x0- xxa0=intervalize $ Bernsteinp (2,1) $ UV.fromList [-x0a,1-x0a] :: Bernsteinp (Int,Int) Interval- yya0=intervalize $ Bernsteinp (1,2) $ UV.fromList [-y0a,1-y0a] :: Bernsteinp (Int,Int) Interval- (Matrix2 aa_ ba_ ca_ da_)=inverse ima- xxa=(scale aa_ xxa0)+(scale ba_ yya0)::Bernsteinp (Int,Int) Interval- yya=(scale ca_ xxa0)+(scale da_ yya0)- - xxb0=intervalize $ Bernsteinp (2,1) $ UV.fromList [-x0b,1-x0b]- yyb0=intervalize $ Bernsteinp (1,2) $ UV.fromList [-y0b,1-y0b]- (Matrix2 ab_ bb_ cb_ db_)=inverse $ intervalize mb- xxb=(scale ab_ xxb0)+(scale bb_ yyb0)- yyb=(scale cb_ xxb0)+(scale db_ yyb0)- - c1=Bernsteinp (1,1) $ UV.singleton 1- - eq0=xxa*xxa+yya*yya-c1- eq1=xxb*xxb+yyb*yyb-c1--inter op@(Bezier{cx=bxa,cy=bya,t0=t0a,t1=t1a}) (Bezier{cx=xb,cy=yb,t0=t0b,t1=t1b})=- - let p0=(promote 1 $ intervalize bxa)-(promote 2 $ intervalize xb) :: Bernsteinp (Int,Int) Interval- p1=(promote 1 $ intervalize bya)-(promote 2 $ intervalize yb) :: Bernsteinp (Int,Int) Interval- thrx=1e-2- solutions _ []=[]- solutions thr boxes@(_:_)=- let sol0=concatMap (solve thr (V.fromList [p0,p1])) boxes - - (correct,toRefine)=partition (\(u,v,_,_)->- let (xu,yu)=evalCurve op (Interval u u)- (xv,yv)=evalCurve op (Interval v v)- in- (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0- in- correct++(solutions (thr/2) toRefine)- in- solutions 1e-2 [(t0a,t1a,t0b,t1b)]---inter cir@(Circle{}) bez@(Bezier{})=map (\(u,v,w,x)->(w,x,u,v)) $ inter bez cir--inter bez@(Bezier{}) cir@(Circle{})=- let xx=(intervalize $ cx bez)-(intervalize $ constant $ cx0 cir)- yy=(intervalize $ cy bez)-(intervalize $ constant $ cy0 cir)- (Matrix2 a b c d)=inverse $ intervalize $ matrix cir- xx0=scale a xx+scale b yy- yy0=scale c xx+scale d yy- - thrx=1e-5- - solutions _ []=[]- solutions thr boxes@(_:_)=- let sol0=concatMap (solve thr (V.singleton (xx0*xx0+yy0*yy0-(constant 1)))) boxes- - (correct,toRefine)=partition (\(u,v)->- let (xu,yu)=evalCurve bez (Interval u u)- (xv,yv)=evalCurve bez (Interval v v)- in- (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0- in- correct++(solutions (thr/2) toRefine)- in- foldl (\l (u,v)->- let (Interval xl xu,Interval yl yu)=evalCurve bez (Interval u v) in- case angle (Interval xl xu) (Interval yl yu) cir of- Just (Interval a0 a1)->- (u,v,a0,a1):l- Nothing->l- ) [] $!- solutions (1e-2) [(t0 bez,t1 bez)]--inter bez@(Bezier{}) off@(Offset{})=map (\(u,v,w,x)->(w,x,u,v)) $ inter off bez---inter off@(Offset{}) bez@(Bezier{})=- - let thr=1e-2- thrx=1e-5- solutions _ []=[]- solutions thr boxes@(_:_)=- let sol0=concatMap (solve thr (V.fromList [eq0,eq1,eq2])) boxes- - (correct,toRefine)=partition (\(u,v,_,_,_,_)->- let (xu,yu)=evalCurve off (Interval u u)- (xv,yv)=evalCurve off (Interval v v)- in- (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0- in- correct++(solutions (thr/2) toRefine)- in- map (\(a,b,c,d,_,_)->(a,b,c,d)) $ solutions 1e-2 $- [(0,1,0,1,0,1)::(Double,Double,Double,Double,Double,Double)]- where- - bxp=intervalize $ cx off- byp=intervalize $ cy off- bxq=intervalize $ cx bez- byq=intervalize $ cy bez- - bxp'=derivate bxp- byp'=derivate byp- bxp3=promote 1 bxp- byp3=promote 1 byp- bxq3=promote 2 bxq- byq3=promote 2 byq- - mp@(Matrix2 ap bp cp dp)=intervalize $ matrix $ off- (Matrix2 ap_ bp_ cp_ dp_)=inverse mp- bXp'=promote 1 $ (scale ap_ bxp')+(scale bp_ byp')- bYp'=promote 1 $ (scale cp_ bxp')+(scale dp_ byp')-- omp@(Bernsteinp _ omegap)=(bXp'*bXp')+(bYp'*bYp')- au=- let au_=sqrt $ UV.minimum $ UV.map ilow omegap in- max 0 $ fpred au_- bu=- let bu_=sqrt $ UV.maximum $ UV.map iup omegap in- fsucc bu_- - alphau=Bernsteinp (1,1,2) $ UV.fromList [Interval au au,Interval bu bu]- eq0=bxp3*alphau + (scale ap bYp') - (scale bp bXp') - bxq3- eq1=byp3*alphau + (scale cp bYp') - (scale dp bXp') - byq3- eq2=alphau*alphau-omp- ---angle::Interval->Interval->Curve->Maybe Interval-angle x y (Circle { cx0,cy0,matrix,t0,t1 })=- let vx=x-interval cx0- vy=y-interval cy0- - Matrix2 a b c d=inverse $ intervalize matrix- -- L'arithmetique d'intervalles fait un peu n'importe quoi- -- quand le vecteur est trop long. On le raccourcit.- alpha=- let co@(Interval col cou)=a*vx+b*vy- Interval sil siu=c*vx+d*vy- co2=- let (col2,cou2)=if col*col<cou*cou then (col*col,cou*cou) else- (cou*cou,col*col)- in- Interval (fpred col2) (fsucc cou2)- si2=- let (sil2,siu2)=if sil*sil<siu*siu then (sil*sil,siu*siu) else- (siu*siu,sil*sil)- in- Interval (fpred sil2) (fsucc siu2)- coco=co/(sqrt (co2+si2))- ac@(Interval acl acu)=acos $ Interval (max (-1) $ ilow coco) (min 1 $ iup coco)- in- if siu<0 then negate ac else- if sil>=0 then ac else- Interval (negate $ min (abs acl) (abs acu))- (max (abs acl) (abs acu))- up ix- | iup ix<t0 =- up $ ix+(2*interval pi)- | otherwise =- down ix- down ix- | ilow ix>t1 =- down $ ix-(2*interval pi)- | iup ix<t0 =- Nothing- | otherwise =- Just ix- in- up alpha---angle _ _ _=error "angle"---- | Pseudo-distance from a point to a curve. Is the result is--- smaller than 1, the point is inside the curve. If it is greater--- than 1, the point is outside. Else we don't know (as usual with--- interval arithmetic).--distance::Interval->Interval->Curve->Interval-distance x0 y0 (Bezier{..})=- distance x0 y0 (Offset{cx,cy,t0,t1,matrix=Matrix2 1 0 0 1})- -distance x0 y0 (Offset{..})=- let (Matrix2 a b c d)=inverse $ intervalize matrix- vx_=intervalize cx-(constant x0)- vy_=intervalize cy-(constant y0)- vx=scale a vx_+scale b vy_- vy=scale c vx_+scale d vy_- - dist=vx*vx+vy*vy- in- foldl (\di (u,v)->let di'=eval dist (Interval u v) in- if iup di<iup di' then di else di') (Interval (1/0) (1/0)) $- (t0,t0):(t1,t1):(solve 1e-5 (V.singleton (derivate dist)) (t0,t1))- - -distance x1 y1 (Circle{..})=- let (Matrix2 a b c d)=inverse $ intervalize matrix- vx_=x1-Interval cx0 cx0- vy_=y1-Interval cy0 cy0- vx=a*vx_+b*vy_- vy=c*vx_+d*vy_- in- vx*vx+vy*vy- --- | Offsets a given Bezier curve with the given pen matrix. The original--- pen is a circle of radius one, the matrix, if inversible, is applied to it.--offset::Matrix2 Double->Curve->[Curve]-offset m (Bezier{cx=x@(Bernsteinp nx bx),cy=y@(Bernsteinp ny by)})=- if nx <=1 && ny <=1 then- [Circle { cx0=UV.head bx,cy0=UV.head by,t0=ilow 0,t1=iup $ 2*pi,matrix=m }]- else- [ c0,c1,c2,c3 ]- - where- im=intervalize m- (Matrix2 a_ b_ c_ d_)=inverse im- - ibx=intervalize x- iby=intervalize y- - lastCoef (Bernsteinp n c)- | n>=1 = UV.last c- | otherwise = 0- firstCoef (Bernsteinp n c)- | n>=1 = UV.head c- | otherwise = 0- - -- Premiere courbe offset- c0=Offset { cx=x, cy=y, t0=0,t1=1,matrix=m }- - -- Demi-cercle 1- - ibx'=derivate ibx- iby'=derivate iby- - -- Calcul du vecteur tangent au bout du premier-- alpha0=- let xx0=lastCoef ibx'- yy0=lastCoef iby'-- xx0'=a_*xx0+b_*yy0- yy0'=c_*xx0+d_*yy0- norm0=sqrt $ xx0'*xx0'+yy0'*yy0'- - xx'=xx0'/norm0- yy'=yy0'/norm0- in- if ilow xx'>=0 then- -(acos yy')- else- if iup xx'<=0 then- acos yy'- else- let Interval u v=acos yy' in- Interval (negate $ max (abs u) (abs v))- (max (abs u) (abs v))- - - alpha0'=alpha0+interval pi- c1=Circle { cx0=lastCoef x,- cy0=lastCoef y,- t0=ilow alpha0,- t1=iup alpha0',- matrix=m }- - - -- Deuxieme courbe offset- c2=Offset { cx=reorient x,- cy=reorient y,- t0=0,t1=1,- matrix=m }- - -- Deuxieme demi-cercle- alpha1=- let xx0=firstCoef ibx'- yy0=firstCoef iby'-- xx0'=a_*xx0+b_*yy0- yy0'=c_*xx0+d_*yy0- norm0=sqrt $ xx0'*xx0'+yy0'*yy0'- - xx'=xx0'/norm0- yy'=yy0'/norm0- in- if ilow xx'>=0 then- -(acos yy')- else- if iup xx'<=0 then- acos yy'- else- let Interval u v=acos yy' in- Interval (negate $ max (abs u) (abs v))- (max (abs u) (abs v))- - alpha1'=alpha1-pi- c3=Circle { cx0=firstCoef x, - cy0=firstCoef y,- t0=ilow alpha1',- t1=iup alpha1,- matrix=m }- --offset _ _=error "offset : undefined yet for other than Bezier"---rnd::Interval->Double-rnd (Interval a b)=(a+b)/2--derivRoots::Double->Curve->([(Double,Double)],[(Double,Double)])-derivRoots thr (Bezier{..})=- (solve thr (V.singleton $ derivate $ intervalize cx) (t0,t1),- solve thr (V.singleton $ derivate $ intervalize cy) (t0,t1))-derivRoots thr (Offset{..})=- let ix=intervalize cx- iy=intervalize cy- x'=derivate ix- y'=derivate iy- m@(Matrix2 a b c d)=intervalize matrix- (Matrix2 a_ b_ c_ d_)=inverse m- - xx'=(scale a_ x')+(scale b_ y')- yy'=(scale c_ x')+(scale d_ y')- xx''=derivate xx'- yy''=derivate yy'- - omega=xx'*xx'+yy'*yy'- au=- let au_=sqrt $ UV.minimum $ UV.map ilow $ coefs omega in- max 0 $ fpred au_- bu=- let bu_=sqrt $ UV.maximum $ UV.map iup $ coefs omega in- fsucc bu_- alphau=Bernsteinp (1,2) $ UV.fromList [Interval au au,Interval bu bu]- - eqx0=(promote 1 yy'')*(alphau^(2::Int))-(promote 1 $ yy'*yy''+xx'*xx'')- eqy0=(promote 1 $ yy'*yy''+xx'*xx'')-(promote 1 xx'')*(alphau^(2::Int))- - eqx=(promote 1 x')*(alphau^(3::Int))+(scale a eqx0)+(scale b eqy0)- :: Bernsteinp (Int,Int) Interval- eqy=(promote 1 y')*(alphau^(3::Int))+(scale c eqx0)+(scale d eqy0)- - eq=(promote 1 omega)-alphau^(2::Int) :: Bernsteinp (Int,Int) Interval- in- (map (\(u,v,_,_)->(u,v)) $ solve thr (V.fromList [eqx,eq])- ((t0,t1,0,1)::(Double,Double,Double,Double)),- map (\(u,v,_,_)->(u,v)) $ solve thr (V.fromList [eqy,eq])- ((t0,t1,0,1)::(Double,Double,Double,Double)))-derivRoots _ (Circle{..})=- let (Matrix2 a b c d)=intervalize matrix- aa=sqrt $ a*a+b*b- cc=sqrt $ c*c+d*d- sx- | ilow a>=0 = acos $ b/aa- | otherwise = negate $ acos $ b/aa- sy- | ilow c>=0 = acos $ d/cc- | otherwise = negate $ acos $ d/cc- Interval ux vx=(-sx)-pi/2- Interval uy vy=(-sy)-pi/2- in- ([(ux,vx)],[(uy,vy)])---- | The leftmost point on a curve-left::Curve->(Double,Double)-left cur=- let (x,y)=- foldl (\m@(xx,_) (s,t)->- let m'@(xx',_)=evalCurve cur $ interval $ (s+t)/2 in- if ilow xx<ilow xx' then m else m') (1/0,1/0) $- (t0 cur,t0 cur):(t1 cur,t1 cur):(fst $ derivRoots 1e-5 cur)- in- (rnd x,rnd y)--- | The bottommost point on a curve-bottom::Curve->(Double,Double)-bottom cur=- let (x,y)=- foldl (\m@(_,yy) (s,t)->- let m'@(_,yy')=evalCurve cur $ interval $ (s+t)/2 in- if ilow yy<ilow yy' then m else m') (1/0,1/0) $- (t0 cur,t0 cur):(t1 cur,t1 cur):(snd $ derivRoots 1e-5 cur)- in- (rnd x,rnd y)--- | The rightmost point on a curve-right::Curve->(Double,Double)-right cur=- let (x,y)=- foldl (\m@(xx,_) (s,t)->- let m'@(xx',_)=evalCurve cur $ interval $ (s+t)/2 in- if iup xx>iup xx' then m else m') (-1/0,-1/0) $- (t0 cur,t0 cur):(t1 cur,t1 cur):(fst $ derivRoots 1e-5 cur)- in- (rnd x,rnd y)--- | The topmost point on a curve-top::Curve->(Double,Double)-top cur=- let (x,y)=- foldl (\m@(_,yy) (s,t)->- let m'@(_,yy')=evalCurve cur $ interval $ (s+t)/2 in- if iup yy>iup yy' then m else m') (-1/0,-1/0) $- (t0 cur,t0 cur):(t1 cur,t1 cur):(snd $ derivRoots 1e-5 cur)- in- (rnd x,rnd y)-\end{code}
+ Graphics/Typography/Geometry.lhs view
@@ -0,0 +1,89 @@+\begin{code}+{-# OPTIONS -XFlexibleInstances -XNamedFieldPuns #-}+-- | This module contains basic tools for geometric types and functions.+module Graphics.Typography.Geometry (Matrix2(..),+ inverse,rotation,+ Geometric(..),+ leftMost,rightMost,topMost,bottomMost)+ where++import Algebra.Polynomials.Numerical++-- | The type of the transformation matrices used in all geometrical applications.+data Matrix2 a=+ -- | The application of @Matrix2 a b c d@ to vector @(x,y)@ should be+ -- @(ax+by,cx+dy)@.+ Matrix2 a a a a deriving (Show, Read, Eq)++-- | Inverses an inversible matrix. If it is not inversible,+-- The behaviour is undefined.+inverse::(Fractional a, Num a)=>Matrix2 a->Matrix2 a+inverse (Matrix2 a b c d)=+ let det=a*d-c*b in+ Matrix2 (d/det) (-b/det) (-c/det) (a/det)++ ++instance Num a=>Num (Matrix2 a) where+ (+) (Matrix2 a b c d) (Matrix2 e f g h)=+ Matrix2 (a+e) (b+f) (c+g) (d+h)+ (*) (Matrix2 a b c d) (Matrix2 e f g h)=+ Matrix2 (a*e+b*g) (a*f+b*h) (c*e+d*g) (c*f+d*h)+ fromInteger a=Matrix2 (fromInteger a) 0 0 (fromInteger a)+ abs=undefined+ signum=undefined++instance Intervalize Matrix2 where+ intervalize (Matrix2 a b c d)=+ Matrix2 (interval a) (interval b) (interval c) (interval d)++ intersects (Matrix2 a b c d) (Matrix2 a' b' c' d')=+ (intersectsd a a') &&+ (intersectsd b b') &&+ (intersectsd c c') &&+ (intersectsd d d')+ +-- | A class for applying geometric applications to objects+class Geometric g where+ translate::Double->Double->g->g+ apply::Matrix2 Double->g->g++-- | The matrix of a rotation+rotation::Floating a=>a->Matrix2 a+rotation theta=+ let ct=cos theta+ st=sin theta+ in+ Matrix2 ct (-st) st ct++instance Geometric g=>Geometric [g] where+ + translate x y cur=map (translate x y) cur+ apply m cur=map (apply m) cur+ ++-- | @'leftMost' a b@ is the leftmost point between @a@ and @b@.+leftMost::(Double,Double)->(Double,Double)->(Double,Double)+leftMost u@(a,_) v@(b,_)+ | a<b = u+ | otherwise = v+-- | @'rightMost' a b@ is the rightmost point between @a@ and @b@.+rightMost::(Double,Double)->(Double,Double)->(Double,Double)+rightMost u@(a,_) v@(b,_)+ | a<b = v+ | otherwise = u+-- | @'bottomMost' a b@ is the lower point between @a@ and @b@.+bottomMost::(Double,Double)->(Double,Double)->(Double,Double)+bottomMost u@(_,a) v@(_,b)+ | a<b = u+ | otherwise = v+-- | @'topMost' a b@ is the upper point between @a@ and @b@.+topMost::(Double,Double)->(Double,Double)->(Double,Double)+topMost u@(_,a) v@(_,b)+ | a<b = v+ | otherwise = u+++++\end{code}
+ Graphics/Typography/Geometry/Approximation.lhs view
@@ -0,0 +1,170 @@+\begin{code}+{-# OPTIONS -XRecordWildCards -XNamedFieldPuns #-}+-- | This module contains the function to approximate a list of curves with+-- degree 3 Bezier curves, using a least squares method.++module Graphics.Typography.Geometry.Approximation(approximate) where++import qualified Data.Vector.Unboxed as UV+import Graphics.Typography.Geometry.Bezier+import Graphics.Typography.Geometry+import Algebra.Polynomials.Bernstein++import Algebra.Polynomials.Numerical+-- import Debug.Trace+rnd::Interval->Double+rnd (Interval a b)=(a+b)/2+++-- | Approximates a list of 'Curves' with a list of degree 3 Bernstein curves.+approximate::[Curve]->[(Bernsteinp Int Double,Bernsteinp Int Double)]+approximate []=[]+approximate l0@(h0:_)= -- traceShow "starting" $+ let approx::Double->Double->[Curve]->[(Bernsteinp Int Double,Bernsteinp Int Double)]+ approx _ _ []=[]+ approx x0 y0 (cc@(Circle {..}):s)= -- traceShow "circle" $+ let theta=abs $ t1-t0 in+ if theta <= pi/2 then+ let x0_=cos $ theta/2+ y0_=sin $ theta/2+ x1_=(4-x0_)/3+ y1_=(1-x0_)*(3-x0_)/(3*y0_)+ + c0=cos $! theta/2+t0+ s0=sin $! theta/2+t0+ + + px0=c0*x0_ - s0*y0_+ py0=s0*x0_ + c0*y0_+ + px1=c0*x1_ - s0*y1_+ py1=s0*x1_ + c0*y1_+ + px2=c0*x1_ + s0*y1_+ py2=s0*x1_ - c0*y1_+ + -- px3=c0*x0_ + s0*y0_+ -- py3=s0*x0_ - c0*y0_+ + x1=cx0+(a*px0+b*py0)+ y1=cy0+(c*px0+d*py0)+ + (Matrix2 a b c d)=matrix+ x=Bernsteinp 4 $ UV.fromList+ [ x0, -- cx0+(a*px3+b*py3),+ cx0+(a*px2+b*py2),+ cx0+(a*px1+b*py1),+ x1]+ + y=Bernsteinp 4 $ UV.fromList+ [ y0, -- cy0+(c*px3+d*py3),+ cy0+(c*px2+d*py2),+ cy0+(c*px1+d*py1),+ y1 ]+ in+ (x,y):(approx x1 y1 s)+ + else+ let t1'=(t1+t0)/2 in+ approx x0 y0 $ (cc { t1=t1' }):(cc { t0=t1' }):s+{-+ approx x0 y0 (h@(Bezier{}):s)=+ -- incorrect !+ (restriction (cx h) (t0 h) (t1 h),+ restriction (cy h) (t0 h) (t1 h)):+ (approx (UV.last $ coefs $ cx h)+ (UV.last $ coefs $ cy h) s)+-}+ -- Ce qui suit est une methode de moindres carres+ approx x0 y0 (off_:s)= -- traceShow ("offset") $+ -- On commence par chercher les points ou la derivee de la norme+ -- de la tangente est maximale. C'est la qu'on va couper s'il y+ -- a un probleme.+ let bx=restriction (cx off_) (t0 off_) (t1 off_)+ by=restriction (cy off_) (t0 off_) (t1 off_)+ off=off_ { cx=bx,cy=by,t0=0,t1=1 }+ ibx=elevate (bounds by-bounds bx) $ intervalize bx+ iby=elevate (bounds bx-bounds by) $ intervalize by+ points=+ let np=10 in+ map (\x->(x/np,x/np)) [0..np]+ -- Ensuite, moindres carres standard, comme dans Hoschek 88.+ + vx0=ibx?1-ibx?0+ vy0=iby?1-iby?0+ vx1=ibx?(bounds ibx-2)-ibx?(bounds ibx-1)+ vy1=iby?(bounds iby-2)-iby?(bounds iby-1)+ + (wx0,wy0)=evalCurve off 0+ (wx1,wy1)=evalCurve off 1++ wx=Bernsteinp 4 $ UV.fromList [wx0,wx0,wx1,wx1] :: Bernsteinp Int Interval+ wy=Bernsteinp 4 $ UV.fromList [wy0,wy0,wy1,wy1] :: Bernsteinp Int Interval++ bern1=Bernsteinp 4 $ UV.fromList [0,1,0,0] :: Bernsteinp Int Interval+ bern2=Bernsteinp 4 $ UV.fromList [0,0,1,0] :: Bernsteinp Int Interval++ sumAll a b c d x1 y1 ((h1,h2):ss)=+ + let h=Interval h1 h2+ (xi,yi)=evalCurve off h+ + b1=eval bern1 h+ b2=eval bern2 h+ + a'=a + (vx0*vx0+vy0*vy0)*b1*b1+ b'=b + (vx0*vx1 + vy0*vy1)*b1*b2+ c'=c + (vx0*vx1 + vy0*vy1)*b1*b2+ d'=d + (vx1*vx1+vy1*vy1)*b2*b2+ + dx=xi-(eval wx h)+ dy=yi-(eval wy h)+ + x1'=x1 + (vx0*dx + vy0*dy)*b1+ y1'=y1 + (vx1*dx + vy1*dy)*b2+ in+ sumAll a' b' c' d' x1' y1' ss+ + sumAll a b c d x1 y1 []=(a,b,c,d,x1,y1)+ + (ra,rb,rc,rd,rx1,ry1)=sumAll 0 0 0 0 0 0 points+ + (Matrix2 a_ b_ c_ d_)=inverse $ Matrix2 ra rb rc rd+ lambda1=a_*rx1+b_*ry1+ lambda2=c_*rx1+d_*ry1+ + -- On a la courbe optimale. Il faut chercher ou on va couper, maintenant+ xap=Bernsteinp 4 $ UV.fromList [wx0,+ wx0+lambda1*vx0,+ wx1+lambda2*vx1,+ wx1]+ yap=Bernsteinp 4 $ UV.fromList [wy0,+ wy0+lambda1*vy0,+ wy1+lambda2*vy1,+ wy1]+ (err,arg)=foldl (\m (h1,h2)->+ let (xi,yi)=evalCurve off (Interval h1 h2)+ xj=eval xap (Interval h1 h2)+ yj=eval yap (Interval h1 h2)+ in+ max m (iup $ abs $ (xi-xj)*(xi-xj)+(yi-yj)*(yi-yj), (h1+h2)/2))+ (0,0) points+ in+ if err<=0.01 then+ (desintervalize xap,desintervalize yap):(approx (rnd wx1) (rnd wy1) s)+ else+ approx x0 y0 $+ (off { cx=restriction (cx off) 0 arg,+ cy=restriction (cy off) 0 arg }):+ (off { cx=restriction (cx off) arg 1,+ cy=restriction (cy off) arg 1 }):s+ + (x0h,y0h)=evalCurve h0 $ Interval (t0 h0) (t0 h0)+ + in+ approx (rnd x0h) (rnd y0h) l0+ +desintervalize::(Bernsteinp a Interval)->(Bernsteinp a Double)+desintervalize b=b { coefs=UV.map rnd $ coefs b}+ +\end{code}
+ Graphics/Typography/Geometry/Bezier.lhs view
@@ -0,0 +1,749 @@+\documentclass{article}+%include lhs2TeX.fmt+\begin{document}+\begin{code}+{-# OPTIONS -XUnboxedTuples -XBangPatterns -XNamedFieldPuns -XRecordWildCards -XMagicHash -cpp #-}+-- | This module contains the basic functions for manipulating Bezier curves. It is heavily+-- based on the book by N. M. Patrikalakis and T. Maekawa, Shape Interrogation for Computer+-- Aided Design and Manufacturing.++module Graphics.Typography.Geometry.Bezier (+ Curve(..),line,bezier3,+ offset,+ inter,+ evalCurve,distance,+ left,bottom,right,top) where++import Algebra.Polynomials.Bernstein+import Algebra.Polynomials.Numerical+import Graphics.Typography.Geometry+import qualified Data.Vector as V+import qualified Data.Vector.Unboxed as UV+import Data.List (partition,sort)++-- | The type for representing all types of curves.+data Curve=+ Bezier { cx::Bernsteinp Int Double,+ cy::Bernsteinp Int Double, + t0::Double,+ t1::Double }+ + | Offset { cx::Bernsteinp Int Double,+ cy::Bernsteinp Int Double,+ t0::Double,+ t1::Double,+ matrix::Matrix2 Double+ }+ | Circle { cx0::Double,+ cy0::Double,+ t0::Double,+ t1::Double,+ matrix::Matrix2 Double+ }+ deriving (Show)+++-- | The basic constructor for lines : a line is a degree 1 Bezier curve+line::Double->Double->Double->Double->Curve+line px py px' py'=Bezier { cx=Bernsteinp 2 $ UV.fromList [px,px'], + cy=Bernsteinp 2 $ UV.fromList [py,py'],+ t0=0,t1=1 }++-- | A shortcut to define degree 3 Bezier curves from points. If the control+-- points are @a,b,c,d@, the function should be called with+-- @'bezier3' xa ya xb yb xc yc xd yd@.+bezier3::Double->Double->Double->Double->Double->Double->Double->Double->Curve+bezier3 px0 py0 px1 py1 px2 py2 px3 py3=+ Bezier { cx=Bernsteinp 4 $ UV.fromList [px0,px1,px2,px3],+ cy=Bernsteinp 4 $ UV.fromList [py0,py1,py2,py3],+ t0=0,t1=1 }++instance Geometric Curve where+ translate x y cur@(Circle{cx0,cy0})= + cur { cx0=cx0+x,cy0=cy0+y }+ translate x y cur=+ cur { cx=(cx cur) { coefs=UV.map (+x) $ coefs $ cx cur},+ cy=(cy cur) { coefs=UV.map (+y) $ coefs $ cy cur} }++ apply m0@(Matrix2 a b c d) cir@(Circle{cx0,cy0,matrix})=+ cir { cx0=a*cx0+b*cy0, cy0=c*cx0+d*cy0, matrix=m0*matrix }+ apply (Matrix2 a b c d) cur=+ cur { cx=(scale a $ cx cur)+(scale b $ cy cur),+ cy=(scale c $ cx cur)+(scale d $ cy cur) }++++-- | Gives the point corresponding to the given value of the parameter+evalCurve::Curve->Interval->(Interval,Interval)+evalCurve (Offset{..}) t=+ let ix=intervalize cx+ iy=intervalize cy+ xt=eval ix t+ yt=eval iy t+ m@(Matrix2 a b c d)=intervalize matrix+ (Matrix2 a_ b_ c_ d_)=inverse m+ xt0'=eval (derivate ix) t+ yt0'=eval (derivate iy) t+ xt'=a_*xt0' + b_*yt0'+ yt'=c_*xt0' + d_*yt0'+ dd=sqrt $ xt'*xt' + yt'*yt'+ in+ (xt+(a*yt'-b*xt')/dd, yt+(c*yt'-d*xt')/dd)++evalCurve (Circle{..}) alpha=+ let xx=cos alpha+ yy=sin alpha+ (Matrix2 a b c d)=intervalize matrix+ in+ (interval cx0+a*xx+b*yy, interval cy0+c*xx+d*yy)+ +evalCurve (Bezier{..}) t=+ let ix=intervalize cx+ iy=intervalize cy+ xx=eval ix t+ yy=eval iy t+ in+ (xx,yy)+ +data Topo=Dehors | SurLaLigne | Dedans deriving Eq++-- | @'inter' c0 c1@ is a list of all possible points of intersection+-- between curves @c0@ and @c1@ : if @(u,v,w,x)@ is returned by 'inter',+-- then curve @c0@ may intersect with @c1@ between parameter values @u@+-- and @v@, which corresponds to parameter values between @w@ and @x@ for+-- @c1@. The implementation guarantees that all actual solutions are found,+-- but possibly false solutions may also be returned.++inter::Curve->Curve->[((Double,Double,Double,Double))]+inter op@(Offset { cx=bxp_,cy=byp_,matrix=mp,t0=t0a,t1=t1a })+ (Offset { cx=bxq_,cy=byq_,matrix=mq,t0=t0b,t1=t1b })=+ + -- Attention : verifier si c'est la meme generatrice+ let thrx=1e-5+ solutions _ []=[]+ solutions thr boxes@(_:_)=+ let sol0=concatMap (solve thr (V.fromList [eq0,eq1,eq2,eq3])) boxes + + (correct,toRefine)=partition (\(u,v,_,_,_,_,_,_)->+ let (xu,yu)=evalCurve op (Interval u u)+ (xv,yv)=evalCurve op (Interval v v)+ in+ (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0+ in+ correct++(solutions (thr/2) toRefine)+ in+ map (\(u,v,w,x,_,_,_,_)->(u,v,w,x)) $ solutions 1e-2 $+ [(t0a,t1a,t0b,t1b,0,1,0,1)::+ (Double,Double,Double,Double,Double,Double,Double,Double)]++ where+ + + imp@(Matrix2 ap bp cp dp)=intervalize mp+ imq@(Matrix2 aq bq cq dq)=intervalize mq+ (Matrix2 ap_ bp_ cp_ dp_)=inverse imp+ (Matrix2 aq_ bq_ cq_ dq_)=inverse imq++ bxp=intervalize bxp_+ byp=intervalize byp_+ bxq=intervalize bxq_+ byq=intervalize byq_+ + bxp4=promote 1 bxp+ byp4=promote 1 byp+ bxq4=promote 2 bxq+ byq4=promote 2 byq+ ++ bxp'=derivate bxp+ byp'=derivate byp+ bxq'=derivate bxq+ byq'=derivate byq++ bXp'=promote 1 $ (scale ap_ bxp')+(scale bp_ byp')+ bYp'=promote 1 $ (scale cp_ bxp')+(scale dp_ byp')++ bXq'=promote 2 $ (scale aq_ bxq')+(scale bq_ byq')+ bYq'=promote 2 $ (scale cq_ bxq')+(scale dq_ byq')++ bomp@(Bernsteinp _ omegap)=(bXp'*bXp')+(bYp'*bYp') :: Bernsteinp (Int,Int,Int,Int) Interval+ bomq@(Bernsteinp _ omegaq)=(bXq'*bXq')+(bYq'*bYq') :: Bernsteinp (Int,Int,Int,Int) Interval+ + au=+ let au_=sqrt $ UV.minimum $ UV.map ilow omegap in+ max 0 $ fpred au_+ bu=+ let bu_=sqrt $ UV.maximum $ UV.map iup omegap in+ fsucc bu_+ av=+ let av_=sqrt $ UV.minimum $ UV.map ilow omegaq in+ max 0 $ fpred av_+ bv=+ let bv_=sqrt $ UV.maximum $ UV.map iup omegaq in+ fsucc bv_+ + alphau=Bernsteinp (1,1,2,1) $ UV.fromList [Interval au au,Interval bu bu]+ alphav=Bernsteinp (1,1,1,2) $ UV.fromList [Interval av av,Interval bv bv]++ eq0=+ ((bxp4*alphau*alphav) + (scale ap $ bYp'*alphav) - (scale bp $ bXp'*alphav)+ -(bxq4*alphau*alphav) - (scale aq $ bYq'*alphau) + (scale bq $ bXq'*alphau))+ eq1=+ ((byp4*alphau*alphav) + (scale cp $ bYp'*alphav) - (scale dp $ bXp'*alphav)+ -(byq4*alphau*alphav) - (scale cq $ bYq'*alphau) + (scale dq $ bXq'*alphau)) + eq2=bomp-(alphau*alphau)+ eq3=bomq-(alphav*alphav)+ + +++inter b@(Circle{}) a@(Offset{})=+ map (\(i,j,k,l)->(k,l,i,j)) $ inter a b++inter o@(Offset { cx=bxp, cy=byp, matrix=mp })+ cir@(Circle{cx0,cy0,matrix=mq})=++ let ix=intervalize bxp+ iy=intervalize byp+ m@(Matrix2 a b c d)=intervalize mp+ (Matrix2 a_ b_ c_ d_)=inverse m+ x'=derivate ix+ y'=derivate iy+ xx'=(scale a_ x')+(scale b_ y')+ yy'=(scale c_ x')+(scale d_ y')+ omega@(Bernsteinp _ omegap)=xx'*xx'+yy'*yy'+ au=+ let au_=sqrt $ UV.minimum $ UV.map ilow omegap in+ max 0 $ fpred au_+ bu=+ let bu_=sqrt $ UV.maximum $ UV.map iup omegap in+ fsucc bu_+ + alphau=Bernsteinp (1,2) $ UV.fromList [Interval au au,Interval bu bu]+ + lambda=(promote 1 omega) - alphau*alphau+ -- Avant multiplication par M_C^-1+ xx0=(promote 1 $ ix-(intervalize $ constant cx0))*alphau+ +(promote 1 $ scale a yy'-scale b xx')+ yy0=(promote 1 $ iy-(intervalize $ constant cy0))*alphau+ +(promote 1 $ scale c yy'-scale d xx') :: Bernsteinp (Int,Int) Interval+ + (Matrix2 ac_ bc_ cc_ dc_)=inverse $ intervalize mq+ xx1=(scale ac_ xx0)+(scale bc_ yy0)+ yy1=(scale cc_ xx0)+(scale dc_ yy0)+ + eqc=xx1*xx1+yy1*yy1-alphau*alphau+ + thrx=1e-5+ + solutions _ []=[]+ solutions thr boxes@(_:_)=+ let sol0=concatMap (solve thr (V.fromList [eqc,lambda])) boxes + + (correct,toRefine)=partition (\(u,v,_,_)->+ let (xu,yu)=evalCurve o (Interval u u)+ (xv,yv)=evalCurve o (Interval v v)+ in+ (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0+ in+ correct++(solutions (thr/2) toRefine)+ ++ -- Removing false positives by computing the distance to the center of+ -- the circle (this is quite fast).+ + removeFalse cl0 (h@(_,v,_,_):h'@(u',_,_,_):s)=+ let u''=(v+u')/2+ (xu,yu)=evalCurve o (Interval u'' u'')+ Interval dl du=distance xu yu cir+ cl1+ | du<1 = Dedans+ | dl>1 = Dehors+ | otherwise = SurLaLigne+ in+ if cl0/=cl1 then h:(removeFalse cl1 (h':s)) else+ removeFalse cl1 (h':s)+ removeFalse _ l=l+ + initCl=+ let (x0,y0)=evalCurve o (Interval (t0 o) (t0 o)) + Interval dl du=distance x0 y0 cir+ in+ if dl>1 then Dehors else if du<1 then Dedans else SurLaLigne+ in+ foldl (\l (u,v,_,_)->+ let (Interval xl xu,Interval yl yu)=evalCurve o (Interval u v) in+ case angle (Interval xl xu) (Interval yl yu) cir of+ Just (Interval a0 a1)->+ (u,v,a0,a1):l+ Nothing->l+ ) [] $ removeFalse initCl $ sort $ solutions 1e-2 [(t0 o,t1 o,0::Double,1::Double)]+ + +inter a@(Circle{cx0=x0a,cy0=y0a,matrix=ma})+ b@(Circle{cx0=x0b,cy0=y0b,matrix=mb})=+ + if (intervalize ma)`intersects`(intervalize mb) && x0a==x0b && y0a==y0b then+ let up ix@(Interval _ x_) tt0 tt1+ | x_<tt0 =+ up (ix+(2*interval pi)) tt0 tt1+ | otherwise = down ix tt0 tt1+ down ix@(Interval x_ x__) tt0 tt1+ | x_>tt1 =+ down (ix-(2*interval pi)) tt0 tt1+ | x__<tt0 =+ Nothing+ | otherwise =+ Just ix+ + alpha=up (interval $ t0 a) (t0 b) (t1 b)+ beta=up (interval $ t0 b) (t0 b) (t1 b)+ in+ + case (alpha,beta) of+ (Just aa,Just ab)->+ case (up aa (t0 a) (t1 a),+ up ab (t0 a) (t1 a)) of+ + (Just ba,Just bb)+ | ilow aa<=iup ab -> [(ilow ba, iup bb,+ ilow aa, iup ab)]+ | otherwise->+ case (up (interval $ t0 b) (t0 a) (t1 a),+ up (interval $ t1 b) (t0 a) (t1 a)) of+ (Just b0,Just b1)->+ [(ilow b0,iup bb,+ t0 b, iup ab),+ (ilow ba,iup b1,+ ilow aa, t1 b)]+ _->[]+ _->[]+ _->[]++ else+ let thr=1e-5+ solutions=solve thr (V.fromList [eq0,eq1]) (fpred u0,fsucc v0,+ fpred w0,fsucc x0)+ in+ foldl (\l (u,v,w,x)->+ let alpha=angle (Interval u v) (Interval w x) a+ beta=angle (Interval u v) (Interval w x) b+ in+ case alpha of+ Just (Interval a0l a0u)->+ case beta of+ Just (Interval b0l b0u)->(a0l,a0u,b0l,b0u):l+ _->l+ _->l+ ) [] solutions + where+ + ima@(Matrix2 am bm cm dm)=intervalize ma+ + maxa=max (iup $ abs am+abs bm) (iup $ abs cm+abs dm)+ (u0,v0,w0,x0)=(x0a-maxa,x0a+maxa,y0a-maxa,y0a+maxa)+ + -- x-x0+ xxa0=intervalize $ Bernsteinp (2,1) $ UV.fromList [-x0a,1-x0a] :: Bernsteinp (Int,Int) Interval+ yya0=intervalize $ Bernsteinp (1,2) $ UV.fromList [-y0a,1-y0a] :: Bernsteinp (Int,Int) Interval+ (Matrix2 aa_ ba_ ca_ da_)=inverse ima+ xxa=(scale aa_ xxa0)+(scale ba_ yya0)::Bernsteinp (Int,Int) Interval+ yya=(scale ca_ xxa0)+(scale da_ yya0)+ + xxb0=intervalize $ Bernsteinp (2,1) $ UV.fromList [-x0b,1-x0b]+ yyb0=intervalize $ Bernsteinp (1,2) $ UV.fromList [-y0b,1-y0b]+ (Matrix2 ab_ bb_ cb_ db_)=inverse $ intervalize mb+ xxb=(scale ab_ xxb0)+(scale bb_ yyb0)+ yyb=(scale cb_ xxb0)+(scale db_ yyb0)+ + c1=Bernsteinp (1,1) $ UV.singleton 1+ + eq0=xxa*xxa+yya*yya-c1+ eq1=xxb*xxb+yyb*yyb-c1++inter op@(Bezier{cx=bxa,cy=bya,t0=t0a,t1=t1a}) (Bezier{cx=xb,cy=yb,t0=t0b,t1=t1b})=+ + let p0=(promote 1 $ intervalize bxa)-(promote 2 $ intervalize xb) :: Bernsteinp (Int,Int) Interval+ p1=(promote 1 $ intervalize bya)-(promote 2 $ intervalize yb) :: Bernsteinp (Int,Int) Interval+ thrx=1e-2+ solutions _ []=[]+ solutions thr boxes@(_:_)=+ let sol0=concatMap (solve thr (V.fromList [p0,p1])) boxes + + (correct,toRefine)=partition (\(u,v,_,_)->+ let (xu,yu)=evalCurve op (Interval u u)+ (xv,yv)=evalCurve op (Interval v v)+ in+ (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0+ in+ correct++(solutions (thr/2) toRefine)+ in+ solutions 1e-2 [(t0a,t1a,t0b,t1b)]+++inter cir@(Circle{}) bez@(Bezier{})=map (\(u,v,w,x)->(w,x,u,v)) $ inter bez cir++inter bez@(Bezier{}) cir@(Circle{})=+ let xx=(intervalize $ cx bez)-(intervalize $ constant $ cx0 cir)+ yy=(intervalize $ cy bez)-(intervalize $ constant $ cy0 cir)+ (Matrix2 a b c d)=inverse $ intervalize $ matrix cir+ xx0=scale a xx+scale b yy+ yy0=scale c xx+scale d yy+ + thrx=1e-5+ + solutions _ []=[]+ solutions thr boxes@(_:_)=+ let sol0=concatMap (solve thr (V.singleton (xx0*xx0+yy0*yy0-(constant 1)))) boxes+ + (correct,toRefine)=partition (\(u,v)->+ let (xu,yu)=evalCurve bez (Interval u u)+ (xv,yv)=evalCurve bez (Interval v v)+ in+ (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0+ in+ correct++(solutions (thr/2) toRefine)+ in+ foldl (\l (u,v)->+ let (Interval xl xu,Interval yl yu)=evalCurve bez (Interval u v) in+ case angle (Interval xl xu) (Interval yl yu) cir of+ Just (Interval a0 a1)->+ (u,v,a0,a1):l+ Nothing->l+ ) [] $!+ solutions (1e-2) [(t0 bez,t1 bez)]++inter bez@(Bezier{}) off@(Offset{})=map (\(u,v,w,x)->(w,x,u,v)) $ inter off bez+++inter off@(Offset{}) bez@(Bezier{})=+ + let thr=1e-2+ thrx=1e-5+ solutions _ []=[]+ solutions thr boxes@(_:_)=+ let sol0=concatMap (solve thr (V.fromList [eq0,eq1,eq2])) boxes+ + (correct,toRefine)=partition (\(u,v,_,_,_,_)->+ let (xu,yu)=evalCurve off (Interval u u)+ (xv,yv)=evalCurve off (Interval v v)+ in+ (iup $ (xu-xv)^(2::Int)+(yu-yv)^(2::Int))<=thrx) sol0+ in+ correct++(solutions (thr/2) toRefine)+ in+ map (\(a,b,c,d,_,_)->(a,b,c,d)) $ solutions 1e-2 $+ [(0,1,0,1,0,1)::(Double,Double,Double,Double,Double,Double)]+ where+ + bxp=intervalize $ cx off+ byp=intervalize $ cy off+ bxq=intervalize $ cx bez+ byq=intervalize $ cy bez+ + bxp'=derivate bxp+ byp'=derivate byp+ bxp3=promote 1 bxp+ byp3=promote 1 byp+ bxq3=promote 2 bxq+ byq3=promote 2 byq+ + mp@(Matrix2 ap bp cp dp)=intervalize $ matrix $ off+ (Matrix2 ap_ bp_ cp_ dp_)=inverse mp+ bXp'=promote 1 $ (scale ap_ bxp')+(scale bp_ byp')+ bYp'=promote 1 $ (scale cp_ bxp')+(scale dp_ byp')++ omp@(Bernsteinp _ omegap)=(bXp'*bXp')+(bYp'*bYp')+ au=+ let au_=sqrt $ UV.minimum $ UV.map ilow omegap in+ max 0 $ fpred au_+ bu=+ let bu_=sqrt $ UV.maximum $ UV.map iup omegap in+ fsucc bu_+ + alphau=Bernsteinp (1,1,2) $ UV.fromList [Interval au au,Interval bu bu]+ eq0=bxp3*alphau + (scale ap bYp') - (scale bp bXp') - bxq3+ eq1=byp3*alphau + (scale cp bYp') - (scale dp bXp') - byq3+ eq2=alphau*alphau-omp+ +++angle::Interval->Interval->Curve->Maybe Interval+angle x y (Circle { cx0,cy0,matrix,t0,t1 })=+ let vx=x-interval cx0+ vy=y-interval cy0+ + Matrix2 a b c d=inverse $ intervalize matrix+ -- L'arithmetique d'intervalles fait un peu n'importe quoi+ -- quand le vecteur est trop long. On le raccourcit.+ alpha=+ let co@(Interval col cou)=a*vx+b*vy+ Interval sil siu=c*vx+d*vy+ co2=+ let (col2,cou2)=if col*col<cou*cou then (col*col,cou*cou) else+ (cou*cou,col*col)+ in+ Interval (fpred col2) (fsucc cou2)+ si2=+ let (sil2,siu2)=if sil*sil<siu*siu then (sil*sil,siu*siu) else+ (siu*siu,sil*sil)+ in+ Interval (fpred sil2) (fsucc siu2)+ coco=co/(sqrt (co2+si2))+ ac@(Interval acl acu)=acos $ Interval (max (-1) $ ilow coco) (min 1 $ iup coco)+ in+ if siu<0 then negate ac else+ if sil>=0 then ac else+ Interval (negate $ min (abs acl) (abs acu))+ (max (abs acl) (abs acu))+ up ix+ | iup ix<t0 =+ up $ ix+(2*interval pi)+ | otherwise =+ down ix+ down ix+ | ilow ix>t1 =+ down $ ix-(2*interval pi)+ | iup ix<t0 =+ Nothing+ | otherwise =+ Just ix+ in+ up alpha+++angle _ _ _=error "angle"++-- | Pseudo-distance from a point to a curve. Is the result is+-- smaller than 1, the point is inside the curve. If it is greater+-- than 1, the point is outside. Else we don't know (as usual with+-- interval arithmetic).++distance::Interval->Interval->Curve->Interval+distance x0 y0 (Bezier{..})=+ distance x0 y0 (Offset{cx,cy,t0,t1,matrix=Matrix2 1 0 0 1})+ +distance x0 y0 (Offset{..})=+ let (Matrix2 a b c d)=inverse $ intervalize matrix+ vx_=intervalize cx-(constant x0)+ vy_=intervalize cy-(constant y0)+ vx=scale a vx_+scale b vy_+ vy=scale c vx_+scale d vy_+ + dist=vx*vx+vy*vy+ in+ foldl (\di (u,v)->let di'=eval dist (Interval u v) in+ if iup di<iup di' then di else di') (Interval (1/0) (1/0)) $+ (t0,t0):(t1,t1):(solve 1e-5 (V.singleton (derivate dist)) (t0,t1))+ + +distance x1 y1 (Circle{..})=+ let (Matrix2 a b c d)=inverse $ intervalize matrix+ vx_=x1-Interval cx0 cx0+ vy_=y1-Interval cy0 cy0+ vx=a*vx_+b*vy_+ vy=c*vx_+d*vy_+ in+ vx*vx+vy*vy+ +-- | Offsets a given Bezier curve with the given pen matrix. The original+-- pen is a circle of radius one, the matrix, if inversible, is applied to it.++offset::Matrix2 Double->Curve->[Curve]+offset m (Bezier{cx=x@(Bernsteinp nx bx),cy=y@(Bernsteinp ny by)})=+ if nx <=1 && ny <=1 then+ [Circle { cx0=UV.head bx,cy0=UV.head by,t0=ilow 0,t1=iup $ 2*pi,matrix=m }]+ else+ [ c0,c1,c2,c3 ]+ + where+ im=intervalize m+ (Matrix2 a_ b_ c_ d_)=inverse im+ + ibx=intervalize x+ iby=intervalize y+ + lastCoef (Bernsteinp n c)+ | n>=1 = UV.last c+ | otherwise = 0+ firstCoef (Bernsteinp n c)+ | n>=1 = UV.head c+ | otherwise = 0+ + -- Premiere courbe offset+ c0=Offset { cx=x, cy=y, t0=0,t1=1,matrix=m }+ + -- Demi-cercle 1+ + ibx'=derivate ibx+ iby'=derivate iby+ + -- Calcul du vecteur tangent au bout du premier++ alpha0=+ let xx0=lastCoef ibx'+ yy0=lastCoef iby'++ xx0'=a_*xx0+b_*yy0+ yy0'=c_*xx0+d_*yy0+ norm0=sqrt $ xx0'*xx0'+yy0'*yy0'+ + xx'=xx0'/norm0+ yy'=yy0'/norm0+ in+ if ilow xx'>=0 then+ -(acos yy')+ else+ if iup xx'<=0 then+ acos yy'+ else+ let Interval u v=acos yy' in+ Interval (negate $ max (abs u) (abs v))+ (max (abs u) (abs v))+ + + alpha0'=alpha0+interval pi+ c1=Circle { cx0=lastCoef x,+ cy0=lastCoef y,+ t0=ilow alpha0,+ t1=iup alpha0',+ matrix=m }+ + + -- Deuxieme courbe offset+ c2=Offset { cx=reorient x,+ cy=reorient y,+ t0=0,t1=1,+ matrix=m }+ + -- Deuxieme demi-cercle+ alpha1=+ let xx0=firstCoef ibx'+ yy0=firstCoef iby'++ xx0'=a_*xx0+b_*yy0+ yy0'=c_*xx0+d_*yy0+ norm0=sqrt $ xx0'*xx0'+yy0'*yy0'+ + xx'=xx0'/norm0+ yy'=yy0'/norm0+ in+ if ilow xx'>=0 then+ -(acos yy')+ else+ if iup xx'<=0 then+ acos yy'+ else+ let Interval u v=acos yy' in+ Interval (negate $ max (abs u) (abs v))+ (max (abs u) (abs v))+ + alpha1'=alpha1-pi+ c3=Circle { cx0=firstCoef x, + cy0=firstCoef y,+ t0=ilow alpha1',+ t1=iup alpha1,+ matrix=m }+ ++offset _ _=error "offset : undefined yet for other than Bezier"+++rnd::Interval->Double+rnd (Interval a b)=(a+b)/2++derivRoots::Double->Curve->([(Double,Double)],[(Double,Double)])+derivRoots thr (Bezier{..})=+ (solve thr (V.singleton $ derivate $ intervalize cx) (t0,t1),+ solve thr (V.singleton $ derivate $ intervalize cy) (t0,t1))+derivRoots thr (Offset{..})=+ let ix=intervalize cx+ iy=intervalize cy+ x'=derivate ix+ y'=derivate iy+ m@(Matrix2 a b c d)=intervalize matrix+ (Matrix2 a_ b_ c_ d_)=inverse m+ + xx'=(scale a_ x')+(scale b_ y')+ yy'=(scale c_ x')+(scale d_ y')+ xx''=derivate xx'+ yy''=derivate yy'+ + omega=xx'*xx'+yy'*yy'+ au=+ let au_=sqrt $ UV.minimum $ UV.map ilow $ coefs omega in+ max 0 $ fpred au_+ bu=+ let bu_=sqrt $ UV.maximum $ UV.map iup $ coefs omega in+ fsucc bu_+ alphau=Bernsteinp (1,2) $ UV.fromList [Interval au au,Interval bu bu]+ + eqx0=(promote 1 yy'')*(alphau^(2::Int))-(promote 1 $ yy'*yy''+xx'*xx'')+ eqy0=(promote 1 $ yy'*yy''+xx'*xx'')-(promote 1 xx'')*(alphau^(2::Int))+ + eqx=(promote 1 x')*(alphau^(3::Int))+(scale a eqx0)+(scale b eqy0)+ :: Bernsteinp (Int,Int) Interval+ eqy=(promote 1 y')*(alphau^(3::Int))+(scale c eqx0)+(scale d eqy0)+ + eq=(promote 1 omega)-alphau^(2::Int) :: Bernsteinp (Int,Int) Interval+ in+ (map (\(u,v,_,_)->(u,v)) $ solve thr (V.fromList [eqx,eq])+ ((t0,t1,0,1)::(Double,Double,Double,Double)),+ map (\(u,v,_,_)->(u,v)) $ solve thr (V.fromList [eqy,eq])+ ((t0,t1,0,1)::(Double,Double,Double,Double)))+derivRoots _ (Circle{..})=+ let (Matrix2 a b c d)=intervalize matrix+ aa=sqrt $ a*a+b*b+ cc=sqrt $ c*c+d*d+ sx+ | ilow a>=0 = acos $ b/aa+ | otherwise = negate $ acos $ b/aa+ sy+ | ilow c>=0 = acos $ d/cc+ | otherwise = negate $ acos $ d/cc+ Interval ux vx=(-sx)-pi/2+ Interval uy vy=(-sy)-pi/2+ in+ ([(ux,vx)],[(uy,vy)])++-- | The leftmost point on a curve+left::Curve->(Double,Double)+left cur=+ let (x,y)=+ foldl (\m@(xx,_) (s,t)->+ let m'@(xx',_)=evalCurve cur $ interval $ (s+t)/2 in+ if ilow xx<ilow xx' then m else m') (1/0,1/0) $+ (t0 cur,t0 cur):(t1 cur,t1 cur):(fst $ derivRoots 1e-5 cur)+ in+ (rnd x,rnd y)+-- | The bottommost point on a curve+bottom::Curve->(Double,Double)+bottom cur=+ let (x,y)=+ foldl (\m@(_,yy) (s,t)->+ let m'@(_,yy')=evalCurve cur $ interval $ (s+t)/2 in+ if ilow yy<ilow yy' then m else m') (1/0,1/0) $+ (t0 cur,t0 cur):(t1 cur,t1 cur):(snd $ derivRoots 1e-5 cur)+ in+ (rnd x,rnd y)+-- | The rightmost point on a curve+right::Curve->(Double,Double)+right cur=+ let (x,y)=+ foldl (\m@(xx,_) (s,t)->+ let m'@(xx',_)=evalCurve cur $ interval $ (s+t)/2 in+ if iup xx>iup xx' then m else m') (-1/0,-1/0) $+ (t0 cur,t0 cur):(t1 cur,t1 cur):(fst $ derivRoots 1e-5 cur)+ in+ (rnd x,rnd y)+-- | The topmost point on a curve+top::Curve->(Double,Double)+top cur=+ let (x,y)=+ foldl (\m@(_,yy) (s,t)->+ let m'@(_,yy')=evalCurve cur $ interval $ (s+t)/2 in+ if iup yy>iup yy' then m else m') (-1/0,-1/0) $+ (t0 cur,t0 cur):(t1 cur,t1 cur):(snd $ derivRoots 1e-5 cur)+ in+ (rnd x,rnd y)+\end{code}
+ Graphics/Typography/Geometry/Outlines.lhs view
@@ -0,0 +1,291 @@+\begin{code}+{-# OPTIONS -XUnboxedTuples -cpp -XRecordWildCards -XNamedFieldPuns -XBangPatterns -XMagicHash -XScopedTypeVariables #-}+-- | This module contains the necessary calls to the other modules of Metafont'+-- to compute the outlines of a given number of pen strokes. The normal way of+-- using it is by calling 'outlines'. One other possible way would be :+--+-- @+-- let curves=cutAll curvesList in+-- remerge $ contour curves $ intersections curves+-- @++module Graphics.Typography.Geometry.Outlines (cutAll, intersections, contour, remerge, outlines) where++import Algebra.Polynomials.Bernstein+import Algebra.Polynomials.Numerical +import Graphics.Typography.Geometry.Bezier+import Graphics.Typography.Geometry+import Data.List (sort)+import qualified Data.Map as M+import qualified Data.Vector as V++import Control.Parallel++(!)::V.Vector a->Int->a+(!)=(V.!)++-- | Cuts a curve into a list of consecutive non-selfintersecting curves.+cutNoSelf::Curve->[Curve]+cutNoSelf c@(Circle{})=[c]+cutNoSelf bez@(Bezier{..})=+ let ix=intervalize cx+ dx=derivate ix+ solutions=+ sort $ filter (\(s,t)->(ilow $ eval ix (Interval s s))*+ (iup $ eval ix (Interval t t)) <= 0) $+ solve 1e-10 (V.singleton dx) (t0,t1)+ roots lastU []=+ if lastU>=t1 then+ []+ else+ [bez { t0=lastU }]+ roots lastU (u:s)+ | u<=lastU = roots lastU s -- on ne coupe pas au debut+ | otherwise =+ (bez { t0=lastU, t1=u }):+ (roots u s)+ in+ roots t0 $ map (\(s,t)->(s+t)/2) solutions+ +cutNoSelf off@(Offset{..})= -- offset+ let thr=1e-2+ ix=intervalize cx+ iy=intervalize cy+ x'=derivate ix+ y'=derivate iy+ (Matrix2 a b c d)=intervalize matrix+ (Matrix2 a_ b_ c_ d_)=inverse $ intervalize matrix+ + xx'=(scale a_ x')+(scale b_ y')+ yy'=(scale c_ x')+(scale d_ y')+ + xx''=derivate xx'+ yy''=derivate yy'+ + evalC (t::Interval)=+ let norm=sqrt $ (eval xx' t)*(eval xx' t)+(eval yy' t)*(eval yy' t)+ derx=(eval yy'' t)/norm - + ((eval yy' t)*((eval xx' t)*(eval xx'' t)++ (eval yy' t)*(eval yy'' t)))/(norm*norm*norm)+ dery=(eval xx'' t)/norm - + ((eval xx' t)*((eval xx' t)*(eval xx'' t)++ (eval yy' t)*(eval yy'' t)))/(norm*norm*norm)+ in+ ((eval x' t)+(a*derx-b*dery), (eval y' t)+(c*derx-d*dery))+ + zerosx=+ let verif t lastxx+ | t>=t1 = []+ | otherwise =+ let (xx,_)=evalC (Interval t t) in+ if (iup $ xx*lastxx)<=0 then+ t:verif (t+thr) xx+ else+ verif (t+thr) xx+ + + + in+ verif t0 $ fst $ evalC (Interval t0 t0)+ + roots lastU []=+ if lastU>=t1 then+ []+ else+ [off { t0=lastU }]+ roots lastU (u:s)+ | u<=lastU = roots lastU s -- on ne coupe pas au debut+ | otherwise =+ (off { t0=lastU, t1=u }):+ (roots u s)+ in+ roots t0 zerosx++-- | @'cutAll' curves@ is the array of all the curves, cut such that+-- each part does not intersect itself.+cutAll::[[Curve]]->V.Vector (V.Vector Curve)+cutAll l=V.fromList $ map (\c->V.fromList $ concatMap cutNoSelf c) l+++data Topology=Dedans | SurLaLigne | Dehors deriving (Eq, Ord, Show)++minsert::Ord a=>a->b->M.Map a [b]->M.Map a [b]+minsert x y m=M.insertWith' (++) x [y] m++munion::Ord a=>M.Map a [b]->M.Map a [b]->M.Map a [b]+munion=M.unionWith (++)+ + +mdeleteFindMin::Ord a=>M.Map a [b]->(Maybe (a,b),M.Map a [b])+mdeleteFindMin m=+ if M.null m then+ (Nothing, m)+ else+ let ((a,b),m')=M.deleteFindMin m in+ case b of+ []->mdeleteFindMin m'+ (h:s)->(Just (a,h), if null s then m' else M.insert a s m')+++-- | Computes the intersections between any pair of curves given+-- as input, in parallel in GHC using @+RTS -N@.+intersections::V.Vector (V.Vector Curve)->+ M.Map (Int,Int,Double) [(Int,Int,Double,Double)]+intersections curves=+ let interAll ci cj+ | ci>=V.length curves = M.empty+ | cj>=V.length curves = interAll (ci+1) (ci+1)+ | otherwise = + -- traceShow (ci,i,cj,j) $+ let next=interAll ci (cj+1)+ inters+ | ci==cj =+ V.ifoldl'+ (\s0 i curvei->+ V.ifoldl' + (\s1 j curvej->+ foldl (\s2 (ti,ti',tj,tj')->+ minsert (ci,i,ti) (cj,j+i+1,tj,tj') $+ minsert (cj,j+i+1,tj) (ci,i,ti,ti') $ s2) s1 $+ inter curvei curvej+ )+ s0 $ V.drop (i+1) (curves!cj)+ ) M.empty $ V.take (V.length (curves!ci)-1) (curves!ci)+ | otherwise = + V.ifoldl'+ (\s0 i curvei->+ V.ifoldl'+ (\s1 j curvej->+ foldl (\s2 (ti,ti',tj,tj')->+ minsert (ci,i,ti) (cj,j,tj,tj') $+ minsert (ci,i,ti') (cj,j,tj,tj') $+ minsert (cj,j,tj) (ci,i,ti,ti') $+ minsert (cj,j,tj') (ci,i,ti,ti') $ s2) s1 $+ inter curvei curvej+ )+ s0 (curves!cj)+ ) M.empty $ V.take (V.length (curves!ci)-1) (curves!ci)+ in+ (next`par`inters)`seq`+ (next`munion`inters)+ in+ interAll 0 0+ +-- | 'contour' takes the curves and the intersections computed as in 'intersections',+-- and outputs a list of all simple closed paths defined by the curves in the input.+contour::V.Vector (V.Vector Curve)->+ M.Map (Int,Int,Double) [(Int,Int,Double,Double)]->+ [[(Int,Int,Double,Double)]]+contour curves inters0=+ + let allPaths inters1 passages1=+ let (first,inters2)=mdeleteFindMin inters1 in+ case first of+ Nothing->[]+ Just ((ci0,i0,ti0),(cj0,j0,tj0a,tj0b))->+ --traceShow ("new path",pi0,pj0) $+ let walk ci i tia tib inters passages=+ --traceShow ("point",ci,i,tia,tib) $ traceShow (inters) $+ let (a,b)=M.split (ci,i,tib) inters+ (next,b')=mdeleteFindMin b+ in+ case next of+ Nothing-> -- traceShow ("echec 1") $+ ([],a,passages)+ Just ((ci',i',ti'),(cj,j,tja,tjb))+ | ci==ci0 && i==i0 && (ci',i',ti')>=(ci,i,ti0)->+ -- fin du chemin+ ([(ci,i,tia,ti0)],a`munion`b',passages)+ + | ci==ci' && i==i' ->+ let isVisible=+ let tt=(tia+ti')/2+ (xi,yi)=evalCurve (curves!ci!i) (Interval tt tt) + in+ V.foldl (\vis cur->+ vis && + iup (distance xi yi $ (cur!0) {t0=0,t1=1})>=1)+ True curves+ in+ if (not isVisible) then+ --traceShow ("invisible",pi') $+ ([],a`munion`b',passages)+ else+ let alreadyPassed=+ let (_,p1)=M.split (ci,i,ti') passages in+ (not $ M.null p1) &&+ (let ((ci_,i_,_),ti'_)=M.findMin p1 in+ ci_==ci && i_==i && ti'_<=ti')+ in+ if alreadyPassed then+ --traceShow ("already passed",pi') $+ ([],a`munion`b',passages)+ else+ --traceShow ("trying",pi') $+ let (nextPath,nextInters,nextPassages)=+ walk cj j tja tjb (a`munion`b') $+ M.insert (ci,i,ti') tia passages+ in+ if null nextPath then+ walk ci i tia tib (a`munion`b') passages+ else+ ((ci,i,tia,ti'):nextPath,+ nextInters,+ M.insert (ci,i,ti') tia nextPassages)+ | otherwise -> --traceShow ("echec 2",ci',i',ti') $+ ([],inters,passages)+ + (path,inters3,passages1')=walk cj0 j0 tj0a tj0b inters2 passages1+ in + if null path then+ --traceShow ("abandon") $+ allPaths inters3 passages1'+ else+ --traceShow ("reussi") $+ path:(allPaths inters3 passages1')+ in+ allPaths inters0 M.empty+ +-- | 'remerge' takes the curves, the output of 'contour', and outputs+-- the list of "remerged" curves, i.e. where the parts free of self-intersections+-- are glued back to each other.+remerge::V.Vector (V.Vector Curve)->[(Int,Int,Double,Double)]->[Curve]+remerge _ []=[]+remerge curves [(ci,i,ti0,ti1)]=[(curves!ci!i) { t0=ti0,t1=ti1 }]+remerge curves (l@((ci,i,ti0,_):s))=+ + let (cj,j,_,tj1)=last s in+ if ci==cj && j+1==i && tj1==ti0 then+ -- dans ce cas, le dernier est colle au premier+ let takeFirsts []=(# [],[] #)+ takeFirsts ((h@(ci',_,_,_)):ss)+ | ci'==ci = + let (# u,v #)=takeFirsts ss in+ (# h:u, v #)+ | otherwise = (# [],h:ss #)+ (# uu,vv #)=takeFirsts l+ in+ remerge_ $ vv++uu+ else+ remerge_ l+ + where+ remerge_ []=[]+ remerge_ [(cj,j,tj0,tj1)]=[(curves!cj!j) { t0=tj0,t1=tj1 }]+ remerge_ ((cj,j,tj0,tj1):(cck@(ck,k,tk0,_)):ss)+ | cj==ck && k==j+1 && tj1==tk0 =+ let (h':s')=remerge_ $ cck:ss in+ (h' { t0=tj0 }) : s'+ + | otherwise = + ((curves!cj!j) { t0=tj0,t1=tj1 }) : (remerge_ $ cck:ss)+ ++-- | Takes a list of curves, potentially offset, and outputs the relevants part+-- of the outlines.+outlines::[[Curve]]->[[Curve]]+outlines curves=+ let curves'=cutAll curves in+ map (remerge curves') $ contour curves' $ intersections curves'++\end{code}
− Graphics/Typography/Outlines.lhs
@@ -1,291 +0,0 @@-\begin{code}-{-# OPTIONS -XUnboxedTuples -cpp -XRecordWildCards -XNamedFieldPuns -XBangPatterns -XMagicHash -XScopedTypeVariables #-}--- | This module contains the necessary calls to the other modules of Metafont'--- to compute the outlines of a given number of pen strokes. The normal way of--- using it is by calling 'outlines'. One other possible way would be :------ @--- let curves=cutAll curvesList in--- remerge $ contour curves $ intersections curves--- @--module Graphics.Typography.Outlines (cutAll, intersections, contour, remerge, outlines) where--import Algebra.Polynomials.Bernstein-import Algebra.Polynomials.Numerical -import Graphics.Typography.Bezier-import Graphics.Typography-import Data.List (sort)-import qualified Data.Map as M-import qualified Data.Vector as V--import Control.Parallel--(!)::V.Vector a->Int->a-(!)=(V.!)---- | Cuts a curve into a list of consecutive non-selfintersecting curves.-cutNoSelf::Curve->[Curve]-cutNoSelf c@(Circle{})=[c]-cutNoSelf bez@(Bezier{..})=- let ix=intervalize cx- dx=derivate ix- solutions=- sort $ filter (\(s,t)->(ilow $ eval ix (Interval s s))*- (iup $ eval ix (Interval t t)) <= 0) $- solve 1e-10 (V.singleton dx) (t0,t1)- roots lastU []=- if lastU>=t1 then- []- else- [bez { t0=lastU }]- roots lastU (u:s)- | u<=lastU = roots lastU s -- on ne coupe pas au debut- | otherwise =- (bez { t0=lastU, t1=u }):- (roots u s)- in- roots t0 $ map (\(s,t)->(s+t)/2) solutions- -cutNoSelf off@(Offset{..})= -- offset- let thr=1e-2- ix=intervalize cx- iy=intervalize cy- x'=derivate ix- y'=derivate iy- (Matrix2 a b c d)=intervalize matrix- (Matrix2 a_ b_ c_ d_)=inverse $ intervalize matrix- - xx'=(scale a_ x')+(scale b_ y')- yy'=(scale c_ x')+(scale d_ y')- - xx''=derivate xx'- yy''=derivate yy'- - evalC (t::Interval)=- let norm=sqrt $ (eval xx' t)*(eval xx' t)+(eval yy' t)*(eval yy' t)- derx=(eval yy'' t)/norm - - ((eval yy' t)*((eval xx' t)*(eval xx'' t)+- (eval yy' t)*(eval yy'' t)))/(norm*norm*norm)- dery=(eval xx'' t)/norm - - ((eval xx' t)*((eval xx' t)*(eval xx'' t)+- (eval yy' t)*(eval yy'' t)))/(norm*norm*norm)- in- ((eval x' t)+(a*derx-b*dery), (eval y' t)+(c*derx-d*dery))- - zerosx=- let verif t lastxx- | t>=t1 = []- | otherwise =- let (xx,_)=evalC (Interval t t) in- if (iup $ xx*lastxx)<=0 then- t:verif (t+thr) xx- else- verif (t+thr) xx- - - - in- verif t0 $ fst $ evalC (Interval t0 t0)- - roots lastU []=- if lastU>=t1 then- []- else- [off { t0=lastU }]- roots lastU (u:s)- | u<=lastU = roots lastU s -- on ne coupe pas au debut- | otherwise =- (off { t0=lastU, t1=u }):- (roots u s)- in- roots t0 zerosx---- | @'cutAll' curves@ is the array of all the curves, cut such that--- each part does not intersect itself.-cutAll::[[Curve]]->V.Vector (V.Vector Curve)-cutAll l=V.fromList $ map (\c->V.fromList $ concatMap cutNoSelf c) l---data Topology=Dedans | SurLaLigne | Dehors deriving (Eq, Ord, Show)--minsert::Ord a=>a->b->M.Map a [b]->M.Map a [b]-minsert x y m=M.insertWith' (++) x [y] m--munion::Ord a=>M.Map a [b]->M.Map a [b]->M.Map a [b]-munion=M.unionWith (++)- - -mdeleteFindMin::Ord a=>M.Map a [b]->(Maybe (a,b),M.Map a [b])-mdeleteFindMin m=- if M.null m then- (Nothing, m)- else- let ((a,b),m')=M.deleteFindMin m in- case b of- []->mdeleteFindMin m'- (h:s)->(Just (a,h), if null s then m' else M.insert a s m')----- | Computes the intersections between any pair of curves given--- as input, in parallel in GHC using @+RTS -N@.-intersections::V.Vector (V.Vector Curve)->- M.Map (Int,Int,Double) [(Int,Int,Double,Double)]-intersections curves=- let interAll ci cj- | ci>=V.length curves = M.empty- | cj>=V.length curves = interAll (ci+1) (ci+1)- | otherwise = - -- traceShow (ci,i,cj,j) $- let next=interAll ci (cj+1)- inters- | ci==cj =- V.ifoldl'- (\s0 i curvei->- V.ifoldl' - (\s1 j curvej->- foldl (\s2 (ti,ti',tj,tj')->- minsert (ci,i,ti) (cj,j+i+1,tj,tj') $- minsert (cj,j+i+1,tj) (ci,i,ti,ti') $ s2) s1 $- inter curvei curvej- )- s0 $ V.drop (i+1) (curves!cj)- ) M.empty $ V.take (V.length (curves!ci)-1) (curves!ci)- | otherwise = - V.ifoldl'- (\s0 i curvei->- V.ifoldl'- (\s1 j curvej->- foldl (\s2 (ti,ti',tj,tj')->- minsert (ci,i,ti) (cj,j,tj,tj') $- minsert (ci,i,ti') (cj,j,tj,tj') $- minsert (cj,j,tj) (ci,i,ti,ti') $- minsert (cj,j,tj') (ci,i,ti,ti') $ s2) s1 $- inter curvei curvej- )- s0 (curves!cj)- ) M.empty $ V.take (V.length (curves!ci)-1) (curves!ci)- in- (next`par`inters)`seq`- (next`munion`inters)- in- interAll 0 0- --- | 'contour' takes the curves and the intersections computed as in 'intersections',--- and outputs a list of all simple closed paths defined by the curves in the input.-contour::V.Vector (V.Vector Curve)->- M.Map (Int,Int,Double) [(Int,Int,Double,Double)]->- [[(Int,Int,Double,Double)]]-contour curves inters0=- - let allPaths inters1 passages1=- let (first,inters2)=mdeleteFindMin inters1 in- case first of- Nothing->[]- Just ((ci0,i0,ti0),(cj0,j0,tj0a,tj0b))->- --traceShow ("new path",pi0,pj0) $- let walk ci i tia tib inters passages=- --traceShow ("point",ci,i,tia,tib) $ traceShow (inters) $- let (a,b)=M.split (ci,i,tib) inters- (next,b')=mdeleteFindMin b- in- case next of- Nothing-> -- traceShow ("echec 1") $- ([],a,passages)- Just ((ci',i',ti'),(cj,j,tja,tjb))- | ci==ci0 && i==i0 && (ci',i',ti')>=(ci,i,ti0)->- -- fin du chemin- ([(ci,i,tia,ti0)],a`munion`b',passages)- - | ci==ci' && i==i' ->- let isVisible=- let tt=(tia+ti')/2- (xi,yi)=evalCurve (curves!ci!i) (Interval tt tt) - in- V.foldl (\vis cur->- vis && - iup (distance xi yi $ (cur!0) {t0=0,t1=1})>=1)- True curves- in- if (not isVisible) then- --traceShow ("invisible",pi') $- ([],a`munion`b',passages)- else- let alreadyPassed=- let (_,p1)=M.split (ci,i,ti') passages in- (not $ M.null p1) &&- (let ((ci_,i_,_),ti'_)=M.findMin p1 in- ci_==ci && i_==i && ti'_<=ti')- in- if alreadyPassed then- --traceShow ("already passed",pi') $- ([],a`munion`b',passages)- else- --traceShow ("trying",pi') $- let (nextPath,nextInters,nextPassages)=- walk cj j tja tjb (a`munion`b') $- M.insert (ci,i,ti') tia passages- in- if null nextPath then- walk ci i tia tib (a`munion`b') passages- else- ((ci,i,tia,ti'):nextPath,- nextInters,- M.insert (ci,i,ti') tia nextPassages)- | otherwise -> --traceShow ("echec 2",ci',i',ti') $- ([],inters,passages)- - (path,inters3,passages1')=walk cj0 j0 tj0a tj0b inters2 passages1- in - if null path then- --traceShow ("abandon") $- allPaths inters3 passages1'- else- --traceShow ("reussi") $- path:(allPaths inters3 passages1')- in- allPaths inters0 M.empty- --- | 'remerge' takes the curves, the output of 'contour', and outputs--- the list of "remerged" curves, i.e. where the parts free of self-intersections--- are glued back to each other.-remerge::V.Vector (V.Vector Curve)->[(Int,Int,Double,Double)]->[Curve]-remerge _ []=[]-remerge curves [(ci,i,ti0,ti1)]=[(curves!ci!i) { t0=ti0,t1=ti1 }]-remerge curves (l@((ci,i,ti0,_):s))=- - let (cj,j,_,tj1)=last s in- if ci==cj && j+1==i && tj1==ti0 then- -- dans ce cas, le dernier est colle au premier- let takeFirsts []=(# [],[] #)- takeFirsts ((h@(ci',_,_,_)):ss)- | ci'==ci = - let (# u,v #)=takeFirsts ss in- (# h:u, v #)- | otherwise = (# [],h:ss #)- (# uu,vv #)=takeFirsts l- in- remerge_ $ vv++uu- else- remerge_ l- - where- remerge_ []=[]- remerge_ [(cj,j,tj0,tj1)]=[(curves!cj!j) { t0=tj0,t1=tj1 }]- remerge_ ((cj,j,tj0,tj1):(cck@(ck,k,tk0,_)):ss)- | cj==ck && k==j+1 && tj1==tk0 =- let (h':s')=remerge_ $ cck:ss in- (h' { t0=tj0 }) : s'- - | otherwise = - ((curves!cj!j) { t0=tj0,t1=tj1 }) : (remerge_ $ cck:ss)- ---- | Takes a list of curves, potentially offset, and outputs the relevants part--- of the outlines.-outlines::[[Curve]]->[[Curve]]-outlines curves=- let curves'=cutAll curves in- map (remerge curves') $ contour curves' $ intersections curves'--\end{code}
typography-geometry.cabal view
@@ -1,5 +1,5 @@ Name: typography-geometry-Version: 1.0+Version: 1.0.0 Synopsis: Drawings for printed text documents Description: Drawings for printed text documents Category: Typography@@ -14,5 +14,6 @@ tag: 1.0 Library Build-Depends: base<5, vector,polynomials-bernstein,containers,parallel- Exposed-modules: Graphics.Typography, Graphics.Typography.Bezier,- Graphics.Typography.Approximation, Graphics.Typography.Outlines+ Exposed-modules: Graphics.Typography.Geometry, Graphics.Typography.Geometry.Bezier,+ Graphics.Typography.Geometry.Approximation,+ Graphics.Typography.Geometry.Outlines