module Main
where
---import Cube
---import Face
---import Grid
---import Misc (flatten)
---import Point
---import RealFunction
---import Tetrahedron
---import ThreeDimensional
-
-
-
---find_point_value :: RealFunction Point
---find_point_value p = poly p
--- where
--- g0 = make_grid 1 trilinear
--- the_cubes = flatten (cubes g0)
--- good_cubes = filter ((flip contains_point) p) the_cubes
--- target_cube = head good_cubes
--- good_tets = filter ((flip contains_point) p) (tetrahedrons target_cube)
--- target_tetrahedron = head good_tets
--- poly = polynomial target_tetrahedron
+import Data.Array.Repa (
+ DIM2,
+ DIM3,
+ Z(..),
+ (:.)(..),
+ slice,
+ reshape,
+ Any(..),
+ All(..)
+ )
+import qualified Data.Array.Repa as R (map)
+
+import Data.Array.Repa.IO.BMP (writeComponentsToBMP)
+import Data.Word
+import System.Environment (getArgs)
+
+import Grid (make_grid, zoom)
+import MRI
+import Values (read_values_3d, write_values_1d)
+
+mri_shape2d :: DIM2
+mri_shape2d = (Z :. 256*2 :. 256*2)
+
+mri_shape3d :: DIM3
+mri_shape3d = (Z :. 256 :. 256 :. 1)
main :: IO ()
main = do
- putStrLn "Hello, World."
- -- print $ find_point_value (0,0,0)
- -- print $ find_point_value (1,0,0)
- -- print $ find_point_value (2,0,0)
- -- print $ find_point_value (0,1,0)
- -- print $ find_point_value (1,1,0)
- -- print $ find_point_value (2,1,0)
- -- print $ find_point_value (0,2,0)
- -- print $ find_point_value (1,2,0)
- -- print $ find_point_value (2,2,0)
- -- print $ find_point_value (0,0,1)
- -- print $ find_point_value (1,0,1)
- -- print $ find_point_value (2,0,1)
- -- print $ find_point_value (0,1,1)
- -- print $ find_point_value (1,1,1)
- -- print $ find_point_value (2,1,1)
- -- print $ find_point_value (0,2,1)
- -- print $ find_point_value (1,2,1)
- -- print $ find_point_value (2,2,1)
- -- print $ find_point_value (0,0,2)
- -- print $ find_point_value (1,0,2)
- -- print $ find_point_value (2,0,2)
- -- print $ find_point_value (0,1,2)
- -- print $ find_point_value (1,1,2)
- -- print $ find_point_value (2,1,2)
- -- print $ find_point_value (0,2,2)
- -- print $ find_point_value (1,2,2)
- -- print $ find_point_value (2,2,2)
+-- args <- getArgs
+-- let color = head args
+-- let in_file = "./data/MRbrain.40." ++ color
+ let out_file = "MRbrain.50.red.out"
+ arr <- read_word16s
+ let arrBrack = bracket_array arr
+ let arrInv = flip_x $ flip_y arrBrack
+ let arrSlice = slice arrInv (Any :. (50 :: Int) :. All :. All)
+ let arrColor = raw_data_to_color arrSlice
+
+ let arrColor' = reshape mri_shape3d arrColor
+ let rdata = red_dbl_data arrColor'
+ let gdata = green_dbl_data arrColor'
+ let bdata = blue_dbl_data arrColor'
+
+-- mridata <- read_values_3d mri_shape in_file
+
+ let gr = make_grid 1 rdata
+ let routput = zoom gr (2,2,1)
+
+ let gg = make_grid 1 gdata
+ let goutput = zoom gg (2,2,1)
+
+ let gb = make_grid 1 bdata
+ let boutput = zoom gb (2,2,1)
+
+ let routput' = R.map double_to_word8 (reshape mri_shape2d routput)
+ let goutput' = R.map double_to_word8 (reshape mri_shape2d goutput)
+ let boutput' = R.map double_to_word8 (reshape mri_shape2d boutput)
+
+-- write_values_1d output out_file
+ writeComponentsToBMP out_file routput' goutput' boutput'
+ where
+ double_to_word8 :: Double -> Word8
+ double_to_word8 x =
+ if x > 255 then
+ 255 :: Word8
+ else
+ if x < 0 then
+ 0 :: Word8
+ else
+ fromIntegral $ truncate x