mapalgebra 0.2.0 → 0.2.1
raw patch · 4 files changed
+31/−288 lines, 4 filesdep −criteriondep ~massiv
Dependencies removed: criterion
Dependency ranges changed: massiv
Files
- CHANGELOG.md +5/−0
- bench/Bench.hs +0/−258
- lib/Geography/MapAlgebra.hs +14/−18
- mapalgebra.cabal +12/−12
CHANGELOG.md view
@@ -1,3 +1,8 @@+## 0.2.1++- Use certain `massiv` functions more canonically. Should improve performance of+ certain Focal Operations.+ ## 0.2.0 - `massiv-0.4` minimum bound.
− bench/Bench.hs
@@ -1,258 +0,0 @@-{-# LANGUAGE DataKinds #-}-{-# LANGUAGE TypeApplications #-}--module Main ( main ) where--import Criterion.Main-import Data.List.NonEmpty (NonEmpty(..))-import Data.Massiv.Array as A hiding (zipWith)-import qualified Data.Vector as V-import qualified Data.Vector.Unboxed as U-import Geography.MapAlgebra-import GHC.TypeLits-import Graphics.ColorSpace-import qualified Numeric.LinearAlgebra as LA-import Prelude hiding (zipWith)-------main :: IO ()-main = do- img <- fileY "data/gray512.tif"- imgF <- fileY "data/gray512.tif"- rgba <- fileRGB "data/512x512.tif"- a512 <- fileRGB "data/512x512.tif"- -- a1024 <- fileRGB "data/1024x1024.tif"- -- a2048 <- fileRGB "data/2048x2048.tif"- -- a4096 <- fileRGB "data/4096x4096.tif"- -- let aHuge :: Raster S u 46500 46500 Word8- -- aHuge = fromFunction S Par (\(r :. c) -> fromIntegral $ r * c)- defaultMain- [- creation- , io- , colouring img- , massivOps img- , localOps rgba a512 img- , hmatrix- , conversions img- , focalOps img imgF- , compositeOps a512 -- a1024 a2048 a4096- , histograms img- ]--creation :: Benchmark-creation = bgroup "Raster Creation"- [ bench "constant 256x256" $ nf (_array . constantB) 5- , bench "constant 512x512" $ nf (_array . constantB') 5- , bench "fromFunction 256x256" $ nf (_array . functionB) (\(r :. c) -> r * c)- , bench "fromFunction 512x512" $ nf (_array . functionB') (\(r :. c) -> r * c)- , bench "fromVector Unboxed Int - 256x256" $ nf (_array . vectorB) uv- , bench "fromVector Boxed Int - 256x256" $ nf (_array . vectorB') bv ]- where uv = U.fromList ([1..65536] :: [Int])- bv = V.fromList ([1..65536] :: [Int])--io :: Benchmark-io = bgroup "IO"- [ bench "fromRGBA 512x512" $ nfIO (_array . _red <$> fileRGB @Word8 @512 @512 "data/512x512.tif")- , bench "fromRGBA (Word8 -> Double) 512x512" $ nfIO (_array . _red <$> fileRGB @Double @512 @512 "data/512x512.tif")- , bench "fromGray Multiband 512x512" $ nfIO (_array <$> fileY @Word8 @512 @512 "data/512x512.tif")- , bench "fromGray Singleband 512x512" $ nfIO (_array <$> gray "data/gray512.tif") ]--colouring :: Raster S p 512 512 Word8 -> Benchmark-colouring img = bgroup "Colouring"- [ bench "classify 512x512" $ nf (_array . strict S . classify invisible cr . lazy) img- , bench "grayscale 512x512" $ nf (_array . strict S . grayscale . lazy) img ]- where cr = greenRed [25, 50, 75, 100, 125, 150, 175, 200, 225, 255]--massivOps :: Raster S p 512 512 Word8 -> Benchmark-massivOps img = bgroup "Massiv Operations"- [ bench "strict S . lazy" $ nf (_array . strict S . lazy) img- , bench "strict U . lazy" $ nf (_array . strict U . lazy) img- , bench "strict P . lazy" $ nf (_array . strict P . lazy) img ]--localOps :: RGBARaster p 512 512 Word8 -> RGBARaster p 512 512 Double -> Raster S p 512 512 Word8 -> Benchmark-localOps (RGBARaster r g b _) (RGBARaster rF gF bF _) img = bgroup "Local Operations"- [ bench "fmap (+ 17) . lazy" $ nf (_array . strict S . fmap (+ 17) . lazy) img- , bench "zipWith (+)" $ nf (_array . strict S . zipWith (+) r) g- , bench "zipWith (/)" $ nf (_array . strict S . zipWith (/) rF) gF- , bench "(+)" $ nf (_array . strict S . (+ lazy r)) (lazy g)- , bench "(/)" $ nf (_array . strict S . (/ lazy rF)) (lazy gF)- -- , bench "(.+)" $ nf (\g' -> computeAs S $ _array r .+ g') (_array g)- -- , bench "(./)" $ nf (\g' -> computeAs S $ _array rF ./ g') (_array gF)- , bench "lmax" $ nf (_array . strict S . lmax img) img- , bench "lmin" $ nf (_array . strict S . lmin img) img- , bench "lmean (Word8)" $ nf (_array . strict S . lmean @Word8 @Double) rs- , bench "lmean (Double)" $ nf (_array . strict S . lmean @Double @Double) rsF- , bench "lvariety" $ nf (_array . strict S . lvariety) rs- , bench "lmajority" $ nf (_array . strict S . lmajority) rs- , bench "lminority" $ nf (_array . strict S . lminority) rs- , bench "lvariance (Word8)" $ nf (fmap (_array . strict S) . lvariance) rs- , bench "lvariance (Double)" $ nf (fmap (_array . strict S) . lvariance) rsF ]- where rs = lazy r :| [lazy g, lazy b]- rsF = lazy rF :| [lazy gF, lazy bF]--hmatrix :: Benchmark-hmatrix = bgroup "HMatrix"- [ bench "linearSolveLS" $ nf (LA.linearSolveLS zing) (LA.col [8,8,8,8,8,8,8,8,8])- , bench "manual - MxM" $ nf (leftPseudo <>) (LA.col [8,8,8,8,8,8,8,8,8])- , bench "manual - MxV" $ nf (leftPseudo LA.#>) (LA.vector [8,8,8,8,8,8,8,8,8]) ]--conversions :: Raster S p 512 512 Word8 -> Benchmark-conversions img = bgroup "Numeric Conversion"- [ bench "Double -> Double via id" $ nf id tau- , bench "Double -> Double via realToFrac" $ nf (realToFrac @Double @Double) tau- , bench "Word -> Double via realToFrac" $ nf (realToFrac @Word8 @Double) 5- , bench "realToFrac on Raster" $ nf (_array . strict S . fmap (realToFrac @Word8 @Double) . lazy) img ]--focalOps :: Raster S p 512 512 Word8 -> Raster S p 512 512 Double -> Benchmark-focalOps img imgF = bgroup "Focal Operations"- [ bgroup "fsum"- [ bench "512" $ nf (_array . strict S . fsum) img- -- , bench "46500" $ nf (_array . strict S . fsum) huge- -- , bench "maybing" $ nf (_array . strict B . fmap Just . lazy) img- -- , bench "nodata" $ nf (_array . strict S . nodatafsum) img- ]- , bgroup "fmean"- [ bench "Word8" $ nf (_array . strict S . fmean @Word8 @Double) img- , bench "Double" $ nf (_array . strict S . fmean @Double @Double) imgF- ]- , bench "fmax" $ nf (_array . strict S . fmax) img- , bench "fmin" $ nf (_array . strict S . fmin) img- , bench "fmajority" $ nf (_array . strict S . fmajority) img- , bench "fminority" $ nf (_array . strict S . fminority) img- , bench "fvariety" $ nf (_array . strict S . fvariety) img- , bench "fpercentage" $ nf (_array . strict S . fpercentage) img- , bench "fpercentile" $ nf (_array . strict S . fpercentile) img- , bench "flinkage" $ nf (_array . strict S . flinkage) img- , bench "flength" $ nf (_array . strict S . flength . flinkage) img- , bench "fpartition" $ nf (_array . strict B . fpartition) img- , bench "fshape" $ nf (_array . strict B . fshape) img- , bench "ffrontage" $ nf (_array . strict S . ffrontage . fshape) img- , bench "farea" $ nf (_array . strict S . farea . fshape) img- , bgroup "fvolume"- [ bench "Word8 -> Double" $ nf (_array . strict S . fvolume . strict S . fmap (realToFrac @Word8 @Double) . lazy) img- , bench "Double" $ nf (_array . strict S . fvolume) imgF- ]- , bgroup "fgradient"- [ bench "Word8" $ nf (_array . strict S . fgradient . strict S . fmap wtod . lazy) img- , bench "Double" $ nf (_array . strict S . fgradient) imgF- ]- , bgroup "faspect"- [ bench "Unsafe (Word8)" $ nf (_array . strict S . faspect' . strict S . fmap wtod . lazy) img- , bench "Unsafe (Double)" $ nf (_array .strict S . faspect') imgF- , bench "Safe (Word8)" $ nf (_array . strict B . faspect . strict S . fmap wtod . lazy) img- , bench "Safe (Double)" $ nf (_array . strict B . faspect) imgF- ]- , bgroup "fdownstream"- [ bench "Word8" $ nf (_array . strict S . fdownstream . strict S . fmap wtod . lazy) img- , bench "Double" $ nf (_array. strict S . fdownstream) imgF- ]- , bgroup "fupstream"- [ bench "Word8" $ nf (_array . strict S . fupstream . strict S . fdownstream . strict S . fmap wtod . lazy) img- , bench "Double" $ nf (_array . strict S . fupstream . strict S . fdownstream) imgF- ]- ]--wtod :: Word8 -> Double-wtod = realToFrac--compositeOps :: RGBARaster p 512 512 Double- -- -> RGBARaster p 1024 1024 Double- -- -> RGBARaster p 2048 2048 Double- -- -> RGBARaster p 4096 4096 Double- -> Benchmark-compositeOps i@(RGBARaster r g _ _) =- bgroup "Composite Operations"- [ bench "NDVI" $ nf (_array . strict S . ndvi (lazy g)) (lazy r)- , bench "EVI (512)" $ nf (_array . strict S . evi) i- -- , bench "EVI (1024)" $ nf (_array . strict S . evi) a1024- -- , bench "EVI (2048)" $ nf (_array . strict S . evi) a2048- -- , bench "EVI (4096)" $ nf (_array . strict S . evi) a4096- , bench "EVI + Colour" $ nf (_array . strict S . classify invisible cr . evi) i- , bench "EVI + Colour + PNG (D)" $ nf (png . classify invisible cr . evi) i- , bench "EVI + Colour + PNG (S)" $ nf (png . strict S . classify invisible cr . evi) i ]- where cr = greenRed $ fmap (10 ^) ([1..10] :: [Int])--fromRight :: Either a b -> b-fromRight (Right b) = b-fromRight _ = error "Was Left"--constantB :: Int -> Raster S p 256 256 Int-constantB = constant S Par--constantB' :: Int -> Raster S p 512 512 Int-constantB' = constant S Par--functionB :: (Ix2 -> Int) -> Raster S p 256 256 Int-functionB = fromFunction S Par--functionB' :: (Ix2 -> Int) -> Raster S p 512 512 Int-functionB' = fromFunction S Par--vectorB :: U.Vector Int -> Raster U p 256 256 Int-vectorB = fromRight . fromVector Par--vectorB' :: V.Vector Int -> Raster B p 256 256 Int-vectorB' = fromRight . fromVector Par---- doubles :: Raster U p 512 512 Double--- doubles = fromRight . fromVector Par $ U.fromList ([1..262144] :: [Double])--zing :: LA.Matrix Double-zing = LA.matrix 3 [ -0.5, -0.5, 1- , -0.5, 0, 1- , -0.5, 0.5, 1- , 0, -0.5, 1- , 0, 0, 1- , 0, 0.5, 1- , 0.5, -0.5, 1- , 0.5, 0, 1- , 0.5, 0.5, 1 ]--gray :: FilePath -> IO (Raster S p 512 512 Word8)-gray fp = do- i <- fromGray fp- case i of- Left err -> error err- Right img -> pure img--fileY :: (Elevator a, KnownNat r, KnownNat c) => FilePath -> IO (Raster S p r c a)-fileY fp = do- i <- fromGray fp- case i of- Left err -> error err- Right img -> pure img-{-# INLINE fileY #-}--fileRGB :: (Elevator a, KnownNat r, KnownNat c) => FilePath -> IO (RGBARaster p r c a)-fileRGB fp = do- i <- fromRGBA fp- case i of- Left err -> error err- Right img -> pure img-{-# INLINE fileRGB #-}---- | See: https://en.wikipedia.org/wiki/Normalized_difference_vegetation_index-ndvi :: (KnownNat r, KnownNat c) => Raster D p r c Double -> Raster D p r c Double -> Raster D p r c Double-ndvi nir red = (nir - red) / (nir + red)-{-# INLINE ndvi #-}---- | See: https://en.wikipedia.org/wiki/Enhanced_vegetation_index-evi :: (Fractional a, Storable a, KnownNat r, KnownNat c) => RGBARaster p r c a -> Raster D p r c a-evi (RGBARaster r g b _) = 2.5 * (numer / denom)- where nir = lazy g -- fudging it.- numer = nir - lazy r- denom = nir + (6 * lazy r) - (7.5 * lazy b) + 1-{-# INLINE evi #-}---- nodatafsum :: Raster S p 1753 1760 Double -> Raster DW p 1753 1760 Double--- nodatafsum = fmap (maybe 0 ((/9) . getSum)) . fmonoid . strict B . fmap check . lazy--- where check 0 = Nothing--- check n = Just $ Sum n--- {-# INLINE nodatafsum #-}--histograms :: Raster S p 512 512 Word8 -> Benchmark-histograms r = bgroup "Histogram Algorithms"- [ bench "Storable Super Mutable" $ nf (_histogram . histogram) r- ]
lib/Geography/MapAlgebra.hs view
@@ -5,7 +5,6 @@ {-# LANGUAGE FlexibleContexts #-} {-# LANGUAGE FlexibleInstances #-} {-# LANGUAGE GeneralizedNewtypeDeriving #-}-{-# LANGUAGE KindSignatures #-} {-# LANGUAGE Rank2Types #-} {-# LANGUAGE ScopedTypeVariables #-} {-# LANGUAGE TypeFamilies #-}@@ -773,40 +772,37 @@ {-# INLINE zipWith #-} -- | Focal Addition.-fsum :: (Num a, Default a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a-fsum (Raster a) = Raster $ mapStencil (Fill 0) (neighbourhoodStencil f) a- where f nw no ne we fo ea sw so se = nw + no + ne + we + fo + ea + sw + so + se+fsum :: (Num a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a+fsum (Raster a) = Raster $ mapStencil (Fill 0) (sumStencil $ Sz2 3 3) a {-# INLINE fsum #-} -- | Focal Product.-fproduct :: (Num a, Default a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a-fproduct (Raster a) = Raster $ mapStencil (Fill 1) (neighbourhoodStencil f) a- where f nw no ne we fo ea sw so se = nw * no * ne * we * fo * ea * sw * so * se+fproduct :: (Num a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a+fproduct (Raster a) = Raster $ mapStencil (Fill 1) (productStencil $ Sz2 3 3) a+{-# INLINE fproduct #-} -- | Focal Monoid - combine all elements of a neighbourhood via their `Monoid` -- instance. In terms of precedence, the neighbourhood focus is the "left-most", -- and all other elements are "added" to it. -- -- This is not mentioned in GaCM, but seems useful nonetheless.-fmonoid :: (Monoid a, Default a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a-fmonoid (Raster a) = Raster $ mapStencil (Fill mempty) (neighbourhoodStencil f) a- where f nw no ne we fo ea sw so se = fo `mappend` nw `mappend` no `mappend` ne `mappend` we `mappend` ea `mappend` sw `mappend` so `mappend` se+fmonoid :: (Monoid a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a+fmonoid (Raster a) = Raster $ mapStencil (Fill mempty) (foldStencil $ Sz2 3 3) a+{-# INLINE fmonoid #-} -- | Focal Mean.-fmean :: (Real a, Fractional b, Default a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c b-fmean = fmap (\n -> realToFrac n / 9) . fsum+fmean :: (Fractional a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a+fmean (Raster a) = Raster $ mapStencil (Fill 0) (avgStencil $ Sz2 3 3) a {-# INLINE fmean #-} -- | Focal Maximum.-fmax :: (Ord a, Default a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a-fmax (Raster a) = Raster $ mapStencil Edge (neighbourhoodStencil f) a- where f nw no ne we fo ea sw so se = P.max nw . P.max no . P.max ne . P.max we . P.max fo . P.max ea . P.max sw $ P.max so se+fmax :: (Ord a, Bounded a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a+fmax (Raster a) = Raster $ mapStencil Edge (maxStencil $ Sz2 3 3) a {-# INLINE fmax #-} -- | Focal Minimum.-fmin :: (Ord a, Default a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a-fmin (Raster a) = Raster $ mapStencil Edge (neighbourhoodStencil f) a- where f nw no ne we fo ea sw so se = P.min nw . P.min no . P.min ne . P.min we . P.min fo . P.min ea . P.min sw $ P.min so se+fmin :: (Ord a, Bounded a, Manifest u Ix2 a) => Raster u p r c a -> Raster DW p r c a+fmin (Raster a) = Raster $ mapStencil Edge (minStencil $ Sz2 3 3) a {-# INLINE fmin #-} -- | Focal Variety - the number of unique values in each neighbourhood.
mapalgebra.cabal view
@@ -1,6 +1,6 @@ cabal-version: 2.2 name: mapalgebra-version: 0.2.0+version: 0.2.1 synopsis: Efficient, polymorphic Map Algebra. description: Efficient, polymorphic Map Algebra.@@ -39,7 +39,7 @@ build-depends: , base >=4.10 && <4.14 , hmatrix >=0.18 && <0.21- , massiv ^>=0.4+ , massiv ^>=0.4.3 , vector >=0.11 && <0.13 library@@ -68,13 +68,13 @@ , tasty-hunit >=0.9 && <0.11 , tasty-quickcheck >=0.8 && <0.11 -benchmark mapalgebra-bench- import: commons- type: exitcode-stdio-1.0- main-is: Bench.hs- hs-source-dirs: bench- ghc-options: -threaded -O2 -with-rtsopts=-N- build-depends:- , criterion >=1.1 && <1.6- , mapalgebra- , massiv-io ^>=0.1+-- benchmark mapalgebra-bench+-- import: commons+-- type: exitcode-stdio-1.0+-- main-is: Bench.hs+-- hs-source-dirs: bench+-- ghc-options: -threaded -O2 -with-rtsopts=-N+-- build-depends:+-- , criterion >=1.1 && <1.6+-- , mapalgebra+-- , massiv-io ^>=0.1