From 715be016934300f596a11e4fc5b8ca2ec42d6c34 Mon Sep 17 00:00:00 2001 From: Michael Orlitzky Date: Mon, 31 Oct 2011 00:29:59 -0400 Subject: [PATCH] Update to repa3 (Ben Lippmeier). --- makefile | 5 ++++- src/Examples.hs | 12 +++++------ src/Grid.hs | 2 +- src/MRI.hs | 54 ++++++++++++++++++++++++++----------------------- src/Main.hs | 27 ++++++++++++++----------- src/Values.hs | 14 ++++++------- 6 files changed, 62 insertions(+), 52 deletions(-) diff --git a/makefile b/makefile index 2bbe576..5ae5bbb 100644 --- a/makefile +++ b/makefile @@ -17,8 +17,12 @@ OPTIMIZATIONS += -funbox-strict-fields OPTIMIZATIONS += -fexcess-precision OPTIMIZATIONS += -fno-spec-constr-count +REPA_INCLUDES := -ivendor/repa-head/repa/dist/build/ +REPA_INCLUDES += -ivendor/repa-head/repa-io/dist/build + GHC_OPTS += $(OPTIMIZATIONS) GHC_OPTS += $(GHC_WARNINGS) +GHC_OPTS += $(REPA_INCLUDES) GHC_OPTS += -odir $(TMPDIR) GHC_OPTS += -hidir $(TMPDIR) GHC_OPTS += --make @@ -26,7 +30,6 @@ GHC_OPTS += -rtsopts GHC_OPTS += -threaded GHC_OPTS += -o bin/${BIN} - .PHONY : test publish_doc doc src_html hlint $(BIN): src/*.hs diff --git a/src/Examples.hs b/src/Examples.hs index 6300ce8..0e805a4 100644 --- a/src/Examples.hs +++ b/src/Examples.hs @@ -33,7 +33,7 @@ n_cube :: Int -> Repa.DIM3 n_cube n = (Repa.Z Repa.:. n Repa.:. n Repa.:. n) trilinear :: Values3D -trilinear = Repa.fromList (n_cube 3) $ +trilinear = Repa.fromListUnboxed (n_cube 3) $ flatten $ transpose_xz trilinear_list @@ -46,7 +46,7 @@ trilinear_zoom_2_list :: [[[Double]]] trilinear_zoom_2_list = [[[1, 3/2, 2, 5/2, 3], [1, 7/4, 5/2, 13/4, 4], [1, 2, 3, 4, 5], [1, 9/4, 7/2, 19/4, 6], [1, 5/2, 4, 11/2, 7]], [[1, 3/2, 2, 5/2, 3], [1, 15/8, 11/4, 29/8, 9/2], [1, 9/4, 7/2, 19/4, 6], [1, 21/8, 17/4, 47/8, 15/2], [1, 3, 5, 7, 9]], [[1, 3/2, 2, 5/2, 3], [1, 2, 3, 4, 5], [1, 5/2, 4, 11/2, 7], [1, 3, 5, 7, 9], [1, 7/2, 6, 17/2, 11]], [[1, 3/2, 2, 5/2, 3], [1, 17/8, 13/4, 35/8, 11/2], [1, 11/4, 9/2, 25/4, 8], [1, 27/8, 23/4, 65/8, 21/2], [1, 4, 7, 10, 13]], [[1, 3/2, 2, 5/2, 3], [1, 9/4, 7/2, 19/4, 6], [1, 3, 5, 7, 9], [1, 15/4, 13/2, 37/4, 12], [1, 9/2, 8, 23/2, 15]]] trilinear_zoom_2 :: Values3D -trilinear_zoom_2 = Repa.fromList (n_cube 6) $ +trilinear_zoom_2 = Repa.fromListUnboxed (n_cube 6) $ flatten $ transpose_xz trilinear_zoom_2_list @@ -60,7 +60,7 @@ trilinear9x9x9_list = [[[1, 1.5, 2, 2.5, 3, 3.5, 4, 4.5, 5], [1, 1.75, 2.5, 3.25 trilinear9x9x9 :: Values3D -trilinear9x9x9 = Repa.fromList (n_cube 9) $ +trilinear9x9x9 = Repa.fromListUnboxed (n_cube 9) $ flatten $ transpose_xz trilinear9x9x9_list @@ -83,7 +83,7 @@ zeros_list = [ [ [ 0, 0, 0 ], -- No need to transpose_xz this one. zeros :: Values3D -zeros = Repa.fromList (n_cube 3) $ flatten zeros_list +zeros = Repa.fromListUnboxed (n_cube 3) $ flatten zeros_list -- | A 3x3x3 array of numbers, starting at (0,0,0) == 0 and counting @@ -102,7 +102,7 @@ naturals_list = [ [ [ 0, 1, 2 ], [ 24, 25, 26 ]]] naturals :: Values3D -naturals = Repa.fromList (n_cube 3) $ +naturals = Repa.fromListUnboxed (n_cube 3) $ flatten $ transpose_xz naturals_list @@ -119,4 +119,4 @@ twenty_vector = (Repa.Z Repa.:. 1 Repa.:. 20 Repa.:. 1) -- | Used in at least one test where we need a 1x20x1 array. naturals_1d :: Values3D -naturals_1d = Repa.fromList twenty_vector (flatten naturals_1d_list) +naturals_1d = Repa.fromListUnboxed twenty_vector (flatten naturals_1d_list) diff --git a/src/Grid.hs b/src/Grid.hs index 59a9387..dc52482 100644 --- a/src/Grid.hs +++ b/src/Grid.hs @@ -137,7 +137,7 @@ zoom :: Values3D -> ScaleFactor -> Values3D zoom v3d scale_factor | xsize == 0 || ysize == 0 || zsize == 0 = empty3d | otherwise = - R.force $ R.unsafeTraverse v3d transExtent f + R.compute $ R.unsafeTraverse v3d transExtent f where (xsize, ysize, zsize) = dims v3d transExtent = zoom_shape scale_factor diff --git a/src/MRI.hs b/src/MRI.hs index c142dcd..f63c729 100644 --- a/src/MRI.hs +++ b/src/MRI.hs @@ -1,3 +1,4 @@ +{-# LANGUAGE FlexibleContexts #-} -- | The MRI module contains functionsd and definitions relevant (and -- specific) to the MRI data files found at, -- @@ -20,9 +21,10 @@ where import Data.Word import Data.Bits import Data.Array.Repa as R +import Data.Array.Repa.Repr.Unboxed as R import Data.Array.Repa.IO.Binary as R -import Data.Array.Repa.IO.ColorRamp as R -import Data.Array.Repa.IO.BMP as R (writeComponentsToBMP) +import Data.Array.Repa.Algorithms.ColorRamp as R +import Data.Array.Repa.IO.BMP as R (writeImageToBMP) import Values @@ -49,20 +51,21 @@ mri_slice3d = (Z :. 1 :. mri_height :. mri_width) -- | RawData is an array of words (16 bits), as contained in the MRI -- data files. -type RawData sh = Array sh Word16 +type RawData sh = Array U sh Word16 -- | A specialization of the 'RawData' type, to three dimensions. type RawData3D = RawData DIM3 type RGB = (Word8, Word8, Word8) -type ColorData sh = Array sh RGB +type ColorData sh = Array U sh RGB +{-# INLINE read_word16s #-} read_word16s :: FilePath -> IO RawData3D read_word16s path = do arr <- R.readArrayFromStorableFile path mri_shape - arr `deepSeqArray` return () - return arr + arr' <- now $ R.copy arr + return arr' bracket :: Double -> Word16 @@ -81,40 +84,46 @@ flip16 xx = shift xx 8 .|. (shift xx (-8) .&. 0x00ff) -swap_bytes :: (Shape sh) => (RawData sh) -> (RawData sh) +{-# INLINE swap_bytes #-} +swap_bytes :: (Shape sh, Repr r Word16) => Array r sh Word16 + -> Array D sh Word16 swap_bytes arr = - R.force $ R.map flip16 arr + R.map flip16 arr -bracket_array :: (Shape sh) => (Values sh) -> (RawData sh) +bracket_array :: Shape sh => Values sh -> Array D sh Word16 bracket_array arr = - R.force $ R.map bracket arr + R.map bracket arr -round_array :: (Shape sh) => (Values sh) -> (RawData sh) +{-# INLINE round_array #-} +round_array :: Shape sh => Values sh -> Array D sh Word16 round_array arr = - R.force $ R.map round arr + R.map round arr -flip_y :: RawData3D -> RawData3D +flip_y :: Repr r Word16 => Array r DIM3 Word16 -> Array D DIM3 Word16 flip_y arr = - R.force $ R.traverse arr id + R.unsafeTraverse arr id (\get (Z :. z :. y :. x) -> get (Z :. z :. (mri_height - 1) - y :. x)) -flip_x :: RawData3D -> RawData3D +flip_x :: Repr r Word16 => Array r DIM3 Word16 -> Array D DIM3 Word16 flip_x arr = - R.force $ R.traverse arr id + R.unsafeTraverse arr id (\get (Z :. z :. y :. x) -> get (Z :. z :. y :. (mri_width - 1) - x)) + +{-# INLINE write_word16s #-} write_word16s :: (Shape sh) => FilePath -> (RawData sh) -> IO () write_word16s = R.writeArrayToStorableFile + values_to_colors :: (Shape sh) => (Values sh) -> (ColorData sh) values_to_colors arr = - R.force $ R.map (truncate_rgb . ramp_it) arr + R.compute $ R.map (truncate_rgb . ramp_it) arr where ramp_it :: Double -> (Double, Double, Double) ramp_it x = @@ -132,19 +141,14 @@ values_to_colors arr = b' = truncate (b * 255) - -z_slice :: Elt a => Int -> Array DIM3 a -> Array DIM2 a +z_slice :: (R.Unbox a, Repr r a) => Int -> Array r DIM3 a -> Array D DIM2 a z_slice n arr = slice arr (Any :. n :. All :. All) - write_values_slice_to_bitmap :: Values2D -> FilePath -> IO () write_values_slice_to_bitmap v3d path = - R.writeComponentsToBMP path routput goutput boutput + R.writeImageToBMP path colors where arr_bracketed = bracket_array v3d - colors = values_to_colors $ R.map fromIntegral arr_bracketed - routput = R.map (\(red, _, _) -> red) colors - goutput = R.map (\(_, green, _) -> green) colors - boutput = R.map (\(_, _, blue) -> blue) colors + colors = values_to_colors $ R.compute $ R.map fromIntegral arr_bracketed diff --git a/src/Main.hs b/src/Main.hs index 80f5e97..cfd1724 100644 --- a/src/Main.hs +++ b/src/Main.hs @@ -32,11 +32,11 @@ main3d = do let zoom_factor = (scale, scale, scale) let out_file = "output.bin" arr <- read_word16s in_file - let arr' = swap_bytes arr - let arrMRI = R.reshape mri_shape arr' - let dbl_data = R.force $ R.map fromIntegral arrMRI - let output = zoom dbl_data zoom_factor - let word16_output = round_array output + let arr' = swap_bytes arr + let arrMRI = R.reshape mri_shape arr' + let dbl_data = R.compute $ R.map fromIntegral arrMRI + let output = zoom dbl_data zoom_factor + let word16_output = R.compute $ round_array output write_word16s out_file word16_output @@ -46,11 +46,14 @@ main2d = do let scale = read s :: Int let zoom_factor = (1, scale, scale) let out_file = "output.bmp" - arr <- read_word16s in_file - let arr' = swap_bytes arr - let arrInv = flip_x $ flip_y arr' - let arrSlice = z_slice 50 arrInv + arr <- read_word16s in_file + let arrSlice = R.computeUnboxed $ z_slice 50 $ flip_x $ flip_y $ swap_bytes arr let arrSlice' = R.reshape mri_slice3d arrSlice - let dbl_data = R.map fromIntegral arrSlice' - let output = zoom dbl_data zoom_factor - write_values_slice_to_bitmap (z_slice 0 output) out_file + + -- If zoom isn't being inlined we need to extract the slice before hand, + -- and convert it to the require formed. + let dbl_data = R.compute $ R.map fromIntegral arrSlice' + let output = zoom dbl_data zoom_factor + let arrSlice0 = R.computeUnboxed $ z_slice 0 output + + write_values_slice_to_bitmap arrSlice0 out_file diff --git a/src/Values.hs b/src/Values.hs index 8e40212..b5dc15a 100644 --- a/src/Values.hs +++ b/src/Values.hs @@ -12,13 +12,13 @@ module Values ( where import Data.Array.Repa ( - Array, + Array, U, Z(..), (:.)(..), DIM2, DIM3, extent, - fromList, + fromListUnboxed, unsafeIndex, ) @@ -29,9 +29,9 @@ import Test.QuickCheck (Arbitrary(..), Gen, choose, vectorOf) import ScaleFactor (ScaleFactor) -type Values sh = Array sh Double -type Values2D = Values DIM2 -type Values3D = Values DIM3 +type Values sh = Array U sh Double +type Values2D = Values DIM2 +type Values3D = Values DIM3 instance Arbitrary Values3D where @@ -42,12 +42,12 @@ instance Arbitrary Values3D where z_dim <- choose (1, 27) elements <- vectorOf (x_dim * y_dim * z_dim) (arbitrary :: Gen Double) let new_shape = (Z :. x_dim :. y_dim :. z_dim) - let three_d = Data.Array.Repa.fromList new_shape elements + let three_d = Data.Array.Repa.fromListUnboxed new_shape elements return three_d empty3d :: Values3D -empty3d = Data.Array.Repa.fromList (Z :. 0 :. 0 :. 0) [] +empty3d = Data.Array.Repa.fromListUnboxed (Z :. 0 :. 0 :. 0) [] dims :: Values3D -> (Int, Int, Int) -- 2.44.2