diff --git a/patch-image.cabal b/patch-image.cabal --- a/patch-image.cabal +++ b/patch-image.cabal @@ -1,5 +1,5 @@ Name: patch-image -Version: 0.3.2.1 +Version: 0.3.3 License: BSD3 License-File: LICENSE Author: Henning Thielemann @@ -41,7 +41,7 @@ README.md Source-Repository this - Tag: 0.3.2.1 + Tag: 0.3.3 Type: darcs Location: http://hub.darcs.net/thielema/patch-image/ @@ -83,11 +83,12 @@ If flag(llvm) Build-Depends: - knead >=0.3 && <0.4, + knead >=0.4 && <0.5, llvm-extra >=0.8 && <0.9, llvm-tf >=3.1 && <3.2, tfp >=1.0 && <1.1, - JuicyPixels >=2.0 && <3.3, + comfort-array >=0.3 && <0.5, + JuicyPixels >=2.0 && <3.4, dsp >=0.2.4 && <0.3, vector >=0.10 && <0.13, pqueue >=1.2 && <1.5, @@ -137,7 +138,7 @@ accelerate-llvm-ptx >=1.1 && <1.2, accelerate-io >=1.0 && <1.1, accelerate >=1.1 && <1.2, - JuicyPixels >=2.0 && <3.3, + JuicyPixels >=2.0 && <3.4, cassava >=0.4.5 && <0.5.2, dsp >=0.2.4 && <0.3, gnuplot >=0.5 && <0.6, @@ -164,7 +165,7 @@ If flag(buildDraft) Build-Depends: - JuicyPixels >=2.0 && <3.3, + JuicyPixels >=2.0 && <3.4, GeomAlgLib >=0.2 && <0.3, utility-ht >=0.0.1 && <0.1, base >=4 && <5 diff --git a/src/Knead.hs b/src/Knead.hs --- a/src/Knead.hs +++ b/src/Knead.hs @@ -9,15 +9,16 @@ import qualified Knead.CArray as KneadCArray import qualified Complex as Komplex import qualified Degree -import MatchImageBorders (arrayCFromKnead, arrayKneadFromC) +import MatchImageBorders (arrayPairFromVec, arrayVecFromPair) import Arithmetic (guardedPairs, maximum0) import LinearAlgebra ( absolutePositionsFromPairDisplacements, fixAtLeastOnePosition, layoutFromPairDisplacements, fixAtLeastOneAnglePosition, ) +import Knead.CArray (liftCArray) import Knead.Shape - (Size, Vec2(Vec2), Dim1, Dim2, Shape2, Index2, Ix2, Factor2, - verticalVal, horizontalVal) + (Size, Vec2(Vec2), Dim1, Dim2, Shape2, Shape2ZB, Index2, Ix2, Factor2, + verticalSize, verticalVal, horizontalVal) import Degree (Degree(Degree), getDegree) import qualified Math.FFT as FFT @@ -27,15 +28,15 @@ import qualified Data.Array.Knead.Simple.Physical as Phys import qualified Data.Array.Knead.Simple.ShapeDependent as ShapeDep import qualified Data.Array.Knead.Simple.Symbolic as Symb -import qualified Data.Array.Knead.Shape.Nested as Shape +import qualified Data.Array.Knead.Shape as Shape import qualified Data.Array.Knead.Expression as Expr import Data.Array.Knead.Simple.Symbolic ((!)) import Data.Array.Knead.Expression (Exp, (==*), (<*), (<=*), (>=*), (||*), (&&*)) -import Data.Array.IArray (amap) -import Data.Array.CArray (CArray) -import Data.Array.MArray (thaw) +import qualified Data.Array.Comfort.Storable.Mutable as MutArray +import qualified Data.Array.Comfort.Storable.Unchecked as ComfortArray +import qualified Data.Array.Comfort.Shape as ComfortShape import qualified LLVM.Extra.ScalarOrVector as SoV import qualified LLVM.Extra.Arithmetic as LLVMArith @@ -88,6 +89,7 @@ type SmallSize = Word32 +type SmallDim = Shape.ZeroBased SmallSize type Plane = Phys.Array Dim2 type SymbPlane = Symb.Array Dim2 @@ -97,7 +99,10 @@ type YUV a = (a,a,a) shape2 :: (Integral i) => i -> i -> Dim2 -shape2 height width = Vec2 (fromIntegral height) (fromIntegral width) +shape2 height width = + Vec2 + (Shape.ZeroBased $ fromIntegral height) + (Shape.ZeroBased $ fromIntegral width) readImage :: Verbosity -> FilePath -> IO ColorImage8 @@ -115,7 +120,7 @@ (Pic.imageHeight pic) (SV.length dat) return $ - Phys.Array + ComfortArray.Array (shape2 (Pic.imageHeight pic) (Pic.imageWidth pic)) (castForeignPtr $ fst $ SV.unsafeToForeignPtr0 dat) _ -> ioError $ userError "unsupported image type" @@ -127,13 +132,14 @@ Phys.Array sh c -> SV.Vector a vectorStorableFrom castArray img = SV.unsafeFromForeignPtr0 - (castArray $ Phys.buffer img) (fromIntegral $ Shape.size $ Phys.shape img) + (castArray $ ComfortArray.buffer img) + (fromIntegral $ ComfortShape.size $ Phys.shape img) imageFromArray :: (Pic.PixelBaseComponent c ~ a, SV.Storable a) => (ForeignPtr b -> ForeignPtr a) -> Phys.Array Dim2 b -> Pic.Image c imageFromArray castArray img = - let Vec2 height width = Phys.shape img + let Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width) = Phys.shape img in Pic.Image { Pic.imageWidth = fromIntegral width, Pic.imageHeight = fromIntegral height, @@ -152,7 +158,8 @@ colorImageExtent :: ColorImage8 -> (Size, Size) colorImageExtent pic = - case Phys.shape pic of Vec2 height width -> (width, height) + case Phys.shape pic of + Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width) -> (width, height) fromInt :: @@ -238,8 +245,8 @@ -atomDim2 :: Shape2 (Atom i) -atomDim2 = Vec2 atom atom +atomDim2 :: Shape2ZB (Atom i) +atomDim2 = Vec2 (Shape.ZeroBased atom) (Shape.ZeroBased atom) atomIx2 :: Index2 (Atom i) atomIx2 = Vec2 atom atom @@ -247,14 +254,16 @@ atomFactor2 :: Factor2 (Atom i) atomFactor2 = Vec2 atom atom -decomposeDim2 :: Exp (Shape2 i) -> Shape2 (Exp i) -decomposeDim2 = Expr.decompose atomDim2 +decomposeDim2 :: Exp (Shape2ZB i) -> Shape2 (Exp i) +decomposeDim2 sh = + case Expr.decompose atomDim2 sh of + Vec2 (Shape.ZeroBased h) (Shape.ZeroBased w) -> Vec2 h w decomposeFactor2 :: Exp (Factor2 i) -> Factor2 (Exp i) decomposeFactor2 = Expr.decompose atomFactor2 -dim2 :: Exp i -> Exp i -> Exp (Shape2 i) -dim2 h w = Expr.compose (Vec2 h w) +dim2 :: Exp i -> Exp i -> Exp (Shape2ZB i) +dim2 h w = Expr.compose (Vec2 (Shape.ZeroBased h) (Shape.ZeroBased w)) ix2 :: Exp i -> Exp i -> Exp (Index2 i) ix2 y x = Expr.compose (Vec2 y x) @@ -509,7 +518,11 @@ tailArr :: (Symb.C array) => array Dim1 a -> array Dim1 a -tailArr = ShapeDep.backpermute (Expr.max 0 . flip Expr.sub 1) (Expr.add 1) +tailArr = + ShapeDep.backpermute + (Expr.modify (Shape.ZeroBased atom) + (fmap (Expr.max 0 . flip Expr.sub 1))) + (Expr.add 1) differentiate :: (Symb.C array, MultiValue.Additive a) => array Dim1 a -> array Dim1 a @@ -546,7 +559,7 @@ (MultiValue.Field a, MultiValue.Real a, MultiValue.RationalConstant a) => SymbPlane a -> SymbPlane a lowpassVert img = - let height = verticalVal $ Symb.shape img + let height = verticalSize $ Symb.shape img in generate (Symb.shape img) $ Expr.modify atomIx2 $ \(Vec2 y x) -> Arith.smooth3 (img ! ix2 (Expr.max 0 (y-1)) x, @@ -603,7 +616,8 @@ (FFTWReal a, MultiValue.Real a, MultiMem.C a, MultiValue.Field a, MultiValue.RationalConstant a) => Dim2 -> IO (Plane a -> Plane a -> IO (Plane a)) -correlatePadded padExtent@(Vec2 height width) = do +correlatePadded + padExtent@(Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width)) = do let sh = Expr.cons padExtent mergePlanes <- RenderP.run $ \a b -> @@ -622,13 +636,6 @@ mulSpecs =<< liftCArray (FFT.dftN [0,1]) =<< mergePlanes a b -liftCArray :: - (SV.Storable a, SV.Storable b) => - (CArray (Int,Int) a -> CArray (Int,Int) b) -> - Plane a -> IO (Plane b) -liftCArray f a = arrayKneadFromC . f <$> arrayCFromKnead a - - prepareOverlapMatching :: IO (Int -> (Degree Float, ColorImage8) -> IO ((Float, Float), Plane Float)) prepareOverlapMatching = do @@ -636,7 +643,7 @@ hp <- highpassMulti rotat <- RenderP.run $ rotate Arith.vecScalar return $ \radius (angle, img) -> - let Vec2 height width = Phys.shape img + let Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width) = Phys.shape img rot = Degree.cis angle ((left, _right), (top, _bottom)) = Arith.boundingBoxOfRotated rot @@ -690,7 +697,8 @@ Exp Float -> Exp Dim2 -> Exp Dim2 -> SymbPlane Float -> SymbPlane (Float, (Size, Size)) -allOverlapsFromCorrelation (Vec2 height width) minOverlapPortion = +allOverlapsFromCorrelation + (Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width)) minOverlapPortion = \sha shb correlated -> let (Vec2 heighta widtha) = decomposeDim2 sha (Vec2 heightb widthb) = decomposeDim2 shb @@ -763,8 +771,12 @@ Symb.map (/ (fromInt xk * fromInt yk)) . Symb.fold1 Expr.add . ShapeDep.backpermute - (Expr.modify atomDim2 $ \(Vec2 height width) -> - (Vec2 (Expr.idiv height yk) (Expr.idiv width xk), Vec2 yk xk)) + (Expr.modify atomDim2 $ + \(Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width)) -> + (Vec2 + (Shape.ZeroBased (Expr.idiv height yk)) + (Shape.ZeroBased (Expr.idiv width xk)), + Vec2 (Shape.ZeroBased yk) (Shape.ZeroBased xk))) (Expr.modify (atomIx2, atomIx2) $ \(Vec2 yi xi, Vec2 yj xj) -> Vec2 (yi*yk+yj) (xi*xk+xj)) @@ -776,9 +788,12 @@ independent from whether minOverlap is zero or not. -} shrinkFactors :: - (Integral a) => Dim2 -> Float -> Shape2 a -> Shape2 a -> Factor2 a -shrinkFactors (Vec2 heightPad widthPad) minOverlapPortion - (Vec2 heighta widtha) (Vec2 heightb widthb) = + (Integral a) => Dim2 -> Float -> Shape2ZB a -> Shape2ZB a -> Factor2 a +shrinkFactors + (Vec2 (Shape.ZeroBased heightPad) (Shape.ZeroBased widthPad)) + minOverlapPortion + (Vec2 (Shape.ZeroBased heighta) (Shape.ZeroBased widtha)) + (Vec2 (Shape.ZeroBased heightb) (Shape.ZeroBased widthb)) = let minOverlap = Arith.minimumOverlapAbsFromPortion minOverlapPortion (min widtha widthb, min heighta heightb) @@ -794,11 +809,14 @@ since we expect that the user chooses an FFT friendly target size. -} shrinkFactorsAlt :: - (Bits a, Integral a) => Float -> Dim2 -> Shape2 a -> Shape2 a -> Factor2 a -shrinkFactorsAlt minOverlapPortion (Vec2 heightPad widthPad) a b = + (Bits a, Integral a) => + Float -> Dim2 -> Shape2ZB a -> Shape2ZB a -> Factor2 a +shrinkFactorsAlt minOverlapPortion + (Vec2 (Shape.ZeroBased heightPad) (Shape.ZeroBased widthPad)) a b = let (widthc,heightc) = Arith.correlationSize minOverlapPortion $ - map (\(Vec2 height width) -> (width, height)) [a,b] + map (\(Vec2 (Shape.ZeroBased height) (Shape.ZeroBased width)) -> + (width, height)) [a,b] in Vec2 (Arith.divUp heightc $ fromIntegral heightPad) (Arith.divUp widthc $ fromIntegral widthPad) @@ -830,10 +848,12 @@ overlappingArea :: (Ord a, Num a) => - Shape2 a -> - Shape2 a -> + Shape2ZB a -> + Shape2ZB a -> (a, a) -> ((a, a), (a, a), (a, a)) -overlappingArea (Vec2 heighta widtha) (Vec2 heightb widthb) (dx, dy) = +overlappingArea + (Vec2 (Shape.ZeroBased heighta) (Shape.ZeroBased widtha)) + (Vec2 (Shape.ZeroBased heightb) (Shape.ZeroBased widthb)) (dx, dy) = let left = max 0 dx top = max 0 dy right = min widtha (widthb + dx) @@ -850,7 +870,9 @@ -} optimalOverlapBigFine :: Dim2 -> IO (Float -> Plane Float -> Plane Float -> IO (Float, (Size, Size))) -optimalOverlapBigFine padExtent@(Vec2 heightPad widthPad) = do +optimalOverlapBigFine + padExtent@(Vec2 (Shape.ZeroBased heightPad) (Shape.ZeroBased widthPad)) = + do overlap <- optimalOverlap padExtent -- optimalOverlap is compiled again here overlapBig <- optimalOverlapBig padExtent @@ -884,13 +906,17 @@ Dim2 -> Dim2 -> Int -> IO (Float -> Maybe Float -> Plane Float -> Plane Float -> IO [(Float, (Size, Size), (Size, Size))]) -optimalOverlapBigMulti padExtent (Vec2 heightStamp widthStamp) numCorrs = do +optimalOverlapBigMulti padExtent + (Vec2 (Shape.ZeroBased heightStamp) (Shape.ZeroBased widthStamp)) + numCorrs = do shrnk <- RenderP.run $ shrink . decomposeFactor2 optOverlap <- optimalOverlap padExtent overDiff <- overlapDifferenceRun clp <- RenderP.run clip - optOverlapFine <- optimalOverlap $ Vec2 (2*heightStamp) (2*widthStamp) + optOverlapFine <- + optimalOverlap $ + Vec2 (Shape.ZeroBased (2*heightStamp)) (Shape.ZeroBased (2*widthStamp)) let overlapFine minimumOverlap a b anchorA@(leftA, topA) anchorB@(leftB, topB) extent@(width,height) = do let addCoarsePos (score, (xm,ym)) = @@ -1127,7 +1153,7 @@ maskedMinimum :: (Shape.C sh, Symb.C array, MultiValue.Real a) => - array (sh, SmallSize) (Bool, a) -> array sh (Bool, a) + array (sh, SmallDim) (Bool, a) -> array sh (Bool, a) maskedMinimum = Symb.fold1 (maybePlus Expr.min) @@ -1191,7 +1217,7 @@ separateDistanceMap :: (Symb.C array, Shape.C sh, MultiValue.C a) => array sh (bool, ((a, a), (a, a))) -> - array (sh, SmallSize) (bool, a) + array (sh, SmallDim) (bool, a) separateDistanceMap array = outerProduct (Expr.modify2 (atom, ((atom, atom), (atom, atom))) atom $ @@ -1200,7 +1226,8 @@ Expr.ifThenElse (expEven $ Expr.idiv sel 2) (uncurry (Expr.ifThenElse (expEven sel)) horiz) (uncurry (Expr.ifThenElse (expEven sel)) vert))) - array (Symb.lift0 $ Symb.id 4) + array (Symb.lift0 $ Symb.id $ + Expr.compose $ Shape.ZeroBased (4 :: Exp SmallSize)) distanceMapBoxRun :: IO (Dim2 -> Geometry Float -> IO (Plane Word8)) @@ -1219,7 +1246,7 @@ (Symb.C array, Shape.C sh, MultiValue.Field a, MultiValue.NativeFloating a ar, MultiValue.Real a, MultiValue.RationalConstant a) => - array SmallSize (Geometry a) -> + array SmallDim (Geometry a) -> array sh (a,a) -> array sh Bool containedAnywhere geoms array = @@ -1237,7 +1264,7 @@ MultiValue.PseudoRing a, MultiValue.Field a, MultiValue.Real a) => Exp Dim2 -> Exp (Geometry a) -> - Symb.Array SmallSize (Geometry a) -> + Symb.Array SmallDim (Geometry a) -> SymbPlane a distanceMapContained sh this others = let distMap = separateDistanceMap $ distanceMapBox sh this @@ -1278,7 +1305,7 @@ (Shape.C sh, Symb.C array, MultiValue.Real a, MultiValue.Algebraic a, MultiValue.IntegerConstant a) => array sh (a,a) -> - array SmallSize (a,a) -> + array SmallDim (a,a) -> array sh a distanceMapPoints a b = Symb.fold1 Expr.min $ @@ -1325,8 +1352,8 @@ MultiValue.NativeFloating a ar) => Exp Dim2 -> Exp (Geometry a) -> - Symb.Array SmallSize (Geometry a) -> - Symb.Array SmallSize (a, a) -> + Symb.Array SmallDim (Geometry a) -> + Symb.Array SmallDim (a, a) -> SymbPlane a distanceMap sh this others points = Symb.zipWith Expr.min @@ -1365,8 +1392,8 @@ Exp a -> Exp Dim2 -> Exp (Geometry a) -> - Symb.Array SmallSize (Geometry a) -> - Symb.Array SmallSize (a, a) -> + Symb.Array SmallDim (Geometry a) -> + Symb.Array SmallDim (a, a) -> SymbPlane a distanceMapGamma gamma sh this others points = Symb.map (pow gamma) $ distanceMap sh this others points @@ -1455,7 +1482,7 @@ pics <- map (mapPicParam (\(State.NoAngleCorrection, pos) -> pos)) <$> processRotation args - let padSize = fromIntegral $ Option.padSize opt + let padSize = Shape.ZeroBased $ fromIntegral $ Option.padSize opt (maybeAllOverlapsShared, optimalOverlapShared) <- case Just $ Vec2 padSize padSize of Just padExtent -> do @@ -1464,6 +1491,7 @@ Nothing -> do let padExtent = uncurry Vec2 $ swap $ + mapPair (Shape.ZeroBased, Shape.ZeroBased) $ Arith.correlationSize (Option.minimumOverlap opt) $ map (colorImageExtent . snd . picColored) pics overlap <- optimalOverlap padExtent @@ -1675,7 +1703,7 @@ when False $ do notice "write fft" let pic0 : pic1 : _ = map snd rotated - size = Vec2 1024 768 + size = Vec2 (Shape.ZeroBased 1024) (Shape.ZeroBased 768) makeByteImage <- RenderP.run $ \k -> imageByteFromFloat . Symb.map (k*) . Symb.fix runPad <- RenderP.run pad @@ -1799,13 +1827,15 @@ picDiffs <- mapM (flip diff avg) rotMovPics getSnd <- RenderP.run $ Symb.map Expr.snd . Symb.fix lp <- lowpassMulti - masks <- map (amap (Bool8.toBool . fst)) <$> mapM arrayCFromKnead picDiffs + let masks = map (ComfortArray.map fst . arrayPairFromVec) picDiffs let smoothRadius = Option.shapeSmooth opt smoothPicDiffs <- - mapM (arrayCFromKnead <=< lp smoothRadius <=< getSnd) picDiffs + mapM (fmap arrayPairFromVec . lp smoothRadius <=< getSnd) picDiffs (locs, pqueue) <- MatchImageBorders.prepareShaping $ zip masks smoothPicDiffs - counts <- thaw . amap (fromIntegral . fst) =<< arrayCFromKnead sumImg + counts <- + MutArray.thaw . ComfortArray.map (fromIntegral . fst) $ + arrayPairFromVec sumImg shapes <- MatchImageBorders.shapeParts counts locs pqueue let names = map FilePath.takeBaseName paths @@ -1813,7 +1843,8 @@ forM_ (Option.outputShapeHard opt) $ \format -> forM_ (zip names shapes) $ \(name,shape) -> writeGrey (Option.quality opt) (printf format name) $ - arrayKneadFromC $ amap (\b -> if b then 255 else 0) shape + arrayVecFromPair $ + ComfortArray.map (\b -> if Bool8.toBool b then 255 else 0) shape emptyPlainCanv <- emptyCanvas addMasked <- addMaskedToCanvas @@ -1821,16 +1852,14 @@ writeImage (Option.quality opt) path =<< foldM (\canvas (shape, rotMovPic) -> - addMasked rotMovPic - (arrayKneadFromC $ amap Bool8.fromBool shape) - canvas) + addMasked rotMovPic (arrayVecFromPair shape) canvas) emptyPlain (zip shapes rotMovPics) forM_ (Option.outputShaped opt) $ \path -> do smoothShapes <- mapM - (lp smoothRadius . arrayKneadFromC . - amap (fromIntegral . fromEnum)) + (lp smoothRadius . arrayVecFromPair . + ComfortArray.map (fromIntegral . fromEnum)) shapes forM_ (Option.outputShape opt) $ \format -> do makeByteImage <- RenderP.run $ imageByteFromFloat . Symb.fix diff --git a/src/Knead/CArray.hs b/src/Knead/CArray.hs --- a/src/Knead/CArray.hs +++ b/src/Knead/CArray.hs @@ -1,5 +1,7 @@ module Knead.CArray where +import Knead.Shape (Vec2(Vec2), Dim2) + import qualified Complex as Komplex import qualified Math.FFT as FFT @@ -8,13 +10,45 @@ import Foreign.Storable (Storable) import Foreign.Storable.Tuple () +import qualified Data.Array.Comfort.Shape as ComfortShape +import Data.Array.Comfort.Storable.Unchecked (Array(Array)) + +import qualified Data.Array.CArray.Base as CArrayPriv import qualified Data.Array.CArray as CArray -import Data.Array.IArray (amap) +import Data.Array.IArray (Ix, amap, bounds, rangeSize) import Data.Array.CArray (CArray) +import Control.Applicative ((<$>)) + import Data.Complex (Complex((:+)), realPart) import Data.Tuple.HT (mapPair) + + +arrayCFromKnead :: Array Dim2 a -> IO (CArray (Int,Int) a) +arrayCFromKnead + (Array + (Vec2 (ComfortShape.ZeroBased height) (ComfortShape.ZeroBased width)) + fptr) = + CArrayPriv.unsafeForeignPtrToCArray fptr + ((0,0), (fromIntegral height - 1, fromIntegral width - 1)) + +arrayKneadFromC :: + (Storable a) => CArray (Int,Int) a -> Array Dim2 a +arrayKneadFromC carray = + case bounds carray of + ((ly,lx), (uy,ux)) -> + Array + (Vec2 + (ComfortShape.ZeroBased $ fromIntegral $ rangeSize (ly,uy)) + (ComfortShape.ZeroBased $ fromIntegral $ rangeSize (lx,ux))) + (snd $ CArrayPriv.toForeignPtr carray) + +liftCArray :: + (Storable a, Storable b) => + (CArray (Int,Int) a -> CArray (Int,Int) b) -> + Array Dim2 a -> IO (Array Dim2 b) +liftCArray f a = arrayKneadFromC . f <$> arrayCFromKnead a pad :: diff --git a/src/Knead/Shape.hs b/src/Knead/Shape.hs --- a/src/Knead/Shape.hs +++ b/src/Knead/Shape.hs @@ -2,9 +2,11 @@ {-# LANGUAGE EmptyDataDecls #-} module Knead.Shape where -import qualified Data.Array.Knead.Shape.Nested as Shape +import qualified Data.Array.Knead.Shape as Shape import qualified Data.Array.Knead.Expression as Expr +import qualified Data.Array.Comfort.Shape as ComfortShape + import qualified LLVM.Extra.Multi.Value.Memory as MultiMem import qualified LLVM.Extra.Multi.Value as MultiValue import qualified LLVM.Extra.Iterator as Iter @@ -22,6 +24,7 @@ import qualified Control.Monad.HT as Monad import Control.Monad (join) +import Data.Tuple.HT (mapSnd) import Data.Int (Int64) @@ -43,8 +46,8 @@ -} type Size = Int64 type Dim0 = () -type Dim1 = Size -type Dim2 = Shape2 Size +type Dim1 = Shape.ZeroBased Size +type Dim2 = Shape2ZB Size type Ix2 = Index2 Size data Vec2 tag i = Vec2 {vertical, horizontal :: i} @@ -54,6 +57,7 @@ data FactorTag type Shape2 = Vec2 ShapeTag +type Shape2ZB n = Shape2 (Shape.ZeroBased n) type Index2 = Vec2 IndexTag type Factor2 = Vec2 FactorTag @@ -137,27 +141,43 @@ zipShape :: MultiValue.T n -> MultiValue.T n -> MultiValue.T (Vec2 tag n) zipShape y x = MultiValue.compose $ Vec2 y x -instance (tag ~ ShapeTag, Shape.C i) => Shape.C (Vec2 tag i) where +instance (tag ~ ShapeTag, ComfortShape.C i) => ComfortShape.C (Vec2 tag i) where + size (Vec2 n m) = ComfortShape.size n * ComfortShape.size m + +instance + (tag ~ ShapeTag, ComfortShape.Indexed i) => + ComfortShape.Indexed (Vec2 tag i) where type Index (Vec2 tag i) = Index2 (Shape.Index i) + indices (Vec2 n m) = map (uncurry Vec2) $ ComfortShape.indices (n,m) + sizeOffset (Vec2 n m) = + mapSnd (. (\(Vec2 i j) -> (i,j))) $ ComfortShape.sizeOffset (n,m) + inBounds (Vec2 n m) (Vec2 i j) = ComfortShape.inBounds (n,m) (i,j) + +instance (tag ~ ShapeTag, Shape.C i) => Shape.C (Vec2 tag i) where intersectCode a b = case (unzipShape a, unzipShape b) of (Vec2 an am, Vec2 bn bm) -> Monad.lift2 zipShape (Shape.intersectCode an bn) (Shape.intersectCode am bm) - sizeCode nm = + size nm = case unzipShape nm of Vec2 n m -> - join $ Monad.lift2 A.mul (Shape.sizeCode n) (Shape.sizeCode m) - size (Vec2 n m) = Shape.size n * Shape.size m - flattenIndexRec nm ij = - case (unzipShape nm, unzipShape ij) of - (Vec2 n m, Vec2 i j) -> do - (ns, il) <- Shape.flattenIndexRec n i - (ms, jl) <- Shape.flattenIndexRec m j - Monad.lift2 (,) - (A.mul ns ms) - (A.add jl =<< A.mul ms il) + join $ Monad.lift2 A.mul (Shape.size n) (Shape.size m) + sizeOffset nm = + case unzipShape nm of + (Vec2 n m) -> do + (ns, iOffset) <- Shape.sizeOffset n + (ms, jOffset) <- Shape.sizeOffset m + sz <- A.mul ns ms + return + (sz, + \ij -> + case unzipShape ij of + (Vec2 i j) -> do + il <- iOffset i + jl <- jOffset j + A.add jl =<< A.mul ms il) loop code nm = case unzipShape nm of Vec2 n m -> @@ -182,3 +202,8 @@ verticalVal, horizontalVal :: (Expr.Value val) => val (Vec2 tag n) -> val n verticalVal = Expr.lift1 (MultiValue.lift1 vertical) horizontalVal = Expr.lift1 (MultiValue.lift1 horizontal) + +verticalSize, horizontalSize :: + (Expr.Value val) => val (Shape2 (Shape.ZeroBased n)) -> val n +verticalSize = Shape.zeroBasedSize . verticalVal +horizontalSize = Shape.zeroBasedSize . horizontalVal diff --git a/src/MatchImageBorders.hs b/src/MatchImageBorders.hs --- a/src/MatchImageBorders.hs +++ b/src/MatchImageBorders.hs @@ -1,3 +1,4 @@ +{-# LANGUAGE TypeFamilies #-} {- | This is an approach for stitching images at narrow bands along lines of small image differences. @@ -10,69 +11,60 @@ -} module MatchImageBorders where -import Knead.Shape (Vec2(Vec2), Dim2) +import Knead.Shape (Vec2(Vec2), Dim2, Size) -import qualified Data.Array.Knead.Simple.Physical as Phys +import qualified Data.Array.Comfort.Storable.Mutable.Unchecked as MutArray +import qualified Data.Array.Comfort.Storable.Unchecked as Array +import qualified Data.Array.Comfort.Shape as ComfortShape +import Data.Array.Comfort.Storable.Mutable (IOArray) +import Data.Array.Comfort.Storable (Array, (!), (//)) import qualified Data.PQueue.Prio.Max as PQ import qualified Data.Set as Set import Data.PQueue.Prio.Max (MaxPQueue) import Data.Set (Set) -import qualified Data.Array.CArray.Base as CArrayPriv -import Data.Array.IOCArray (IOCArray) -import Data.Array.MArray (readArray, writeArray, freeze, thaw) -import Data.Array.CArray (CArray) -import Data.Array.IArray - (Ix, amap, bounds, range, rangeSize, inRange, (!), (//)) - -import Foreign.Storable (Storable) - +import qualified Data.Bool8 as Bool8 import Data.Traversable (forM) -import Data.Tuple.HT (mapSnd) +import Data.Foldable (forM_) +import Data.Tuple.HT (mapPair, mapSnd) import Data.Maybe (mapMaybe, listToMaybe) import Data.Word (Word8) +import Data.Bool8 (Bool8) import Control.Monad (filterM) import Control.Applicative ((<$>)) -arrayCFromKnead :: Phys.Array Dim2 a -> IO (CArray (Int,Int) a) -arrayCFromKnead (Phys.Array (Vec2 height width) fptr) = - CArrayPriv.unsafeForeignPtrToCArray fptr - ((0,0), (fromIntegral height - 1, fromIntegral width - 1)) +type Z2 i j = (ComfortShape.ZeroBased i, ComfortShape.ZeroBased j) -arrayKneadFromC :: - (Storable a) => CArray (Int,Int) a -> Phys.Array Dim2 a -arrayKneadFromC carray = - case bounds carray of - ((ly,lx), (uy,ux)) -> - Phys.Array - (Vec2 - (fromIntegral (rangeSize (ly,uy))) - (fromIntegral (rangeSize (lx,ux)))) - (snd $ CArrayPriv.toForeignPtr carray) +arrayPairFromVec :: Array Dim2 a -> Array (Z2 Size Size) a +arrayPairFromVec = Array.mapShape (\(Vec2 height width) -> (height, width)) +arrayVecFromPair :: Array (Z2 Size Size) a -> Array Dim2 a +arrayVecFromPair = Array.mapShape (\(height, width) -> (Vec2 height width)) -findBorder :: (Ix i, Enum i, Ix j, Enum j) => CArray (i,j) Bool -> Set (i,j) + +findBorder :: (Integral i, Integral j) => Array (Z2 i j) Bool8 -> Set (i,j) findBorder mask = - let ((yl,xl), (yu,xu)) = bounds mask - revRange (l,u) = [u, pred u .. l] - findLeft y = - listToMaybe $ dropWhile (\x -> not $ mask!(y,x)) $ range (xl,xu) - findRight y = - listToMaybe $ dropWhile (\x -> not $ mask!(y,x)) $ revRange (xl,xu) - findTop x = - listToMaybe $ dropWhile (\y -> not $ mask!(y,x)) $ range (yl,yu) - findBottom x = - listToMaybe $ dropWhile (\y -> not $ mask!(y,x)) $ revRange (yl,yu) + let (yrng, xrng) = Array.shape mask + range sh = ComfortShape.indices sh + revRange sh@(ComfortShape.ZeroBased n) = + take (ComfortShape.size sh) $ drop 1 $ iterate (subtract 1) n + first p = listToMaybe . dropWhile (not . Bool8.toBool . p) + findLeft y = first (\x -> mask!(y,x)) $ range xrng + findRight y = first (\x -> mask!(y,x)) $ revRange xrng + findTop x = first (\y -> mask!(y,x)) $ range yrng + findBottom x = first (\y -> mask!(y,x)) $ revRange yrng in Set.fromList $ - mapMaybe (\y -> (,) y <$> findLeft y) (range (yl,yu)) ++ - mapMaybe (\y -> (,) y <$> findRight y) (range (yl,yu)) ++ - mapMaybe (\x -> flip (,) x <$> findTop x) (range (xl,xu)) ++ - mapMaybe (\x -> flip (,) x <$> findBottom x) (range (xl,xu)) + mapMaybe (\y -> (,) y <$> findLeft y) (range yrng) ++ + mapMaybe (\y -> (,) y <$> findRight y) (range yrng) ++ + mapMaybe (\x -> flip (,) x <$> findTop x) (range xrng) ++ + mapMaybe (\x -> flip (,) x <$> findBottom x) (range xrng) -pqueueFromBorder :: (Ix ix) => CArray ix Float -> Set ix -> MaxPQueue Float ix +pqueueFromBorder :: + (ComfortShape.Indexed sh, ComfortShape.Index sh ~ ix) => + Array sh Float -> Set ix -> MaxPQueue Float ix pqueueFromBorder weights = PQ.fromList . map (\pos -> (weights!pos, pos)) . Set.toList @@ -84,55 +76,62 @@ locBorder = 1 locInside = 2 -prepareLocations :: (Ix ix) => CArray ix Bool -> Set ix -> CArray ix Location +prepareLocations :: + (ComfortShape.Indexed sh, ComfortShape.Index sh ~ ix) => + Array sh Bool8 -> Set ix -> Array sh Location prepareLocations mask border = - amap (\b -> if b then locInside else locOutside) mask + Array.map (\b -> if Bool8.toBool b then locInside else locOutside) mask // map (flip (,) locBorder) (Set.toList border) type Queue i j = - MaxPQueue Float ((IOCArray (i,j) Location, CArray (i,j) Float), (i,j)) + MaxPQueue Float + ((IOArray (Z2 i j) Location, Array (Z2 i j) Float), (i,j)) prepareShaping :: - (Ix i, Enum i, Ix j, Enum j) => - [(CArray (i,j) Bool, CArray (i,j) Float)] -> - IO ([IOCArray (i,j) Location], Queue i j) + (Integral i, Integral j) => + [(Array (Z2 i j) Bool8, Array (Z2 i j) Float)] -> + IO ([IOArray (Z2 i j) Location], Queue i j) prepareShaping maskWeightss = fmap (mapSnd PQ.unions . unzip) $ forM maskWeightss $ \(mask, weights) -> do let border = findBorder mask - locations <- thaw $ prepareLocations mask border + locations <- MutArray.thaw $ prepareLocations mask border return (locations, fmap ((,) (locations, weights)) $ pqueueFromBorder weights border) -shapeParts :: - IOCArray (Int, Int) Int -> - [IOCArray (Int, Int) Location] -> - Queue Int Int -> IO [CArray (Int, Int) Bool] -shapeParts count masks = - let loop :: Queue Int Int -> IO [CArray (Int, Int) Bool] - loop queue = + +loopQueue :: + (Monad m, Ord k) => (a -> m (MaxPQueue k a)) -> MaxPQueue k a -> m () +loopQueue f = + let loop queue = case PQ.maxView queue of - Nothing -> mapM (fmap (amap (/=locOutside)) . freeze) masks - Just (((locs, diffs), pos@(y,x)), remQueue) -> do - n <- readArray count pos - if n<=1 - then loop remQueue - else do - writeArray count pos (n-1) - writeArray locs pos locOutside - envPoss <- - filterM (fmap (locInside ==) . readArray locs) $ - filter (inRange (bounds diffs)) $ - map - (\(dy,dx) -> (y+dy, x+dx)) - [(0,1), (1,0), (0,-1), (-1,0)] - mapM_ - (\envPos -> writeArray locs envPos locBorder) - envPoss - loop $ PQ.union remQueue $ PQ.fromList $ - map (\envPos -> (diffs!envPos, ((locs, diffs), envPos))) envPoss - in loop + Nothing -> return () + Just (first, remQueue) -> loop . PQ.union remQueue =<< f first + in loop + +shapeParts :: + IOArray (Z2 Size Size) Int -> + [IOArray (Z2 Size Size) Location] -> + Queue Size Size -> IO [Array (Z2 Size Size) Bool8] +shapeParts count masks queue = do + flip loopQueue queue $ \((locs, diffs), pos@(y,x)) -> do + n <- MutArray.read count pos + if n<=1 + then return PQ.empty + else do + MutArray.write count pos (n-1) + MutArray.write locs pos locOutside + envPoss <- + filterM (fmap (locInside ==) . MutArray.read locs) $ + filter (ComfortShape.inBounds (Array.shape diffs)) $ + map (mapPair ((y+), (x+))) [(0,1), (1,0), (0,-1), (-1,0)] + forM_ envPoss $ \envPos -> MutArray.write locs envPos locBorder + return $ PQ.fromList $ + map (\envPos -> (diffs!envPos, ((locs, diffs), envPos))) envPoss + + forM masks $ + fmap (Array.map (Bool8.fromBool . (/=locOutside))) . MutArray.freeze