-
Notifications
You must be signed in to change notification settings - Fork 19
Deconvolution algorithm #26
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: master
Are you sure you want to change the base?
Changes from all commits
9c1d31b
08c298f
9b016b3
c2cdcb5
e6b0015
ea50a9b
1509901
4733cf2
2d1de1a
1915739
e1f5626
8757017
aeefdc3
eded16d
4ee7104
5436509
bf8c508
d92fdc4
df4db10
a1f23dd
8fd18c6
258b6ad
fe345d2
22a7003
26c9a06
File filter
Filter by extension
Conversations
Jump to
Diff view
Diff view
There are no files selected for viewing
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,85 @@ | ||
| {-# LANGUAGE RankNTypes #-} | ||
| {-# LANGUAGE FlexibleContexts #-} | ||
| {-# LANGUAGE ScopedTypeVariables #-} | ||
| {-# LANGUAGE BangPatterns #-} | ||
|
|
||
| -- | Adaptive histogram equalization is used to improve contrast in images. | ||
| -- It adjusts image intensity in small regions (neighborhood) in the image. | ||
| module Graphics.Image.Processing.Ahe where | ||
|
|
||
| import Control.Monad (forM_, when) | ||
| import Control.Monad.ST | ||
| import Data.STRef | ||
| import Debug.Trace (trace) | ||
|
|
||
| import Prelude as P hiding (subtract) | ||
| import Graphics.Image.Processing.Filter | ||
| import Graphics.Image.Interface as I | ||
| import Graphics.Image | ||
| import Graphics.Image.Types as IP | ||
| import Graphics.Image.ColorSpace (X) | ||
|
|
||
| -- | Supplementary function for applying border resolution and a general mask. | ||
| simpleFilter :: (Array arr cs e, Array arr X e) => Direction -> Border (Pixel cs e) -> Filter arr cs e | ||
| simpleFilter dir !border = | ||
| Filter (correlate border kernel) | ||
| where | ||
| !kernel = | ||
| case dir of | ||
| Vertical -> fromLists $ [ [ 0, -1, 0 ], [ -1, 4, -1 ], [ 0, -1, 0 ] ] | ||
| Horizontal -> fromLists $ [ [ 0, -1, 0 ], [ -1, 4, -1 ], [ 0, -1, 0 ] ] | ||
|
|
||
| -- | 'ahe' operates on small 'contextual' regions of the image. It enhances the contrast of each | ||
| -- region and this technique works well when the distribution of pixel values is similar throughout | ||
| -- the image. | ||
| -- | ||
| -- The idea is to perform contrast enhancement in 'neighborhood region' of each pixel and the size | ||
| -- of the region is a parameter of the method. It constitutes a characteristic length scale: contrast | ||
| -- at smaller scales is enhanced, while contrast at larger scales is reduced (For general purposes, a size | ||
| -- factor of 5 tends to give pretty good results). | ||
| -- | ||
| -- <<images/yield.jpg>> <<images/yield_ahe.png>> | ||
| -- | ||
| -- Usage : | ||
| -- | ||
| -- >>> img <- readImageY VU "images/yield.jpg" | ||
| -- >>> input1 <- getLine | ||
| -- >>> input2 <- getLine | ||
| -- >>> let thetaSz = (P.read input1 :: Int) | ||
| -- >>> let distSz = (P.read input2 :: Int) | ||
| -- >>> let neighborhoodFactor = (P.read input2 :: Int) | ||
| -- >>> let aheImage :: Image VU RGB Double | ||
| -- >>> aheImage = ahe img thetaSz distSz neighborhoodFactor | ||
| -- >>> writeImage "images/yield_ahe.png" (toImageRGB aheImage) | ||
| -- | ||
| ahe | ||
| :: forall arr e cs . ( MArray arr Y Double, IP.Array arr Y Double, IP.Array arr Y Word16, MArray arr Y Word16, Array arr X Double) | ||
| => Image arr Y Double | ||
| -> Int -- ^ width of output image | ||
| -> Int -- ^ height of output image | ||
| -> Int -- ^ neighborhood size factor | ||
| -> Image arr Y Word16 | ||
| ahe image thetaSz distSz neighborhoodFactor = I.map (fmap toWord16) accBin | ||
| where | ||
| ip = applyFilter (simpleFilter Horizontal Edge) image -- Pre-processing (Border resolution) | ||
| widthMax, var1, heightMax, var2 :: Int | ||
| var1 = ((rows ip) - 1) | ||
| widthMax = ((rows ip) - 1) | ||
| var2 = ((cols ip) - 1) | ||
| heightMax = ((cols ip) - 1) | ||
|
|
||
| accBin :: Image arr Y Word16 | ||
| accBin = runST $ -- Core part of the Algo begins here. | ||
| do arr <- I.new (thetaSz, distSz) -- Create a mutable image with the given dimensions. | ||
| forM_ [0 .. var1] $ \x -> do | ||
| forM_ [0 .. var2] $ \y -> do | ||
| rankRef <- newSTRef (0 :: Int) | ||
| let neighborhood a maxValue = filter (\a -> a >= 0 && a < maxValue) [a-5 .. a+5] | ||
| forM_ (neighborhood x var1) $ \i -> do | ||
| forM_ (neighborhood y var2) $ \j -> do | ||
| when (I.index ip (x, y) > I.index ip (i, j)) $ modifySTRef' rankRef (+1) | ||
| rank <- readSTRef rankRef | ||
| let px = ((rank * 255)) | ||
| I.write arr (x, y) (PixelY (fromIntegral px)) | ||
| freeze arr | ||
|
|
| Original file line number | Diff line number | Diff line change |
|---|---|---|
| @@ -0,0 +1,49 @@ | ||
| {-# LANGUAGE RankNTypes #-} | ||
| {-# LANGUAGE FlexibleContexts #-} | ||
| {-# LANGUAGE ScopedTypeVariables #-} | ||
| {-# LANGUAGE BangPatterns #-} | ||
|
|
||
| module Graphics.Image.Processing.Deconvolution where | ||
|
|
||
| import Graphics.Image as I | ||
| import Graphics.Image.Interface as I | ||
| import Prelude as P | ||
|
|
||
| -- | Supplementary function needed to blur an image and test the algorithm. | ||
| randomBlur :: Image VS X Double | ||
| randomBlur = randomBlur' / scalar (I.sum randomBlur') | ||
| where randomBlur' = fromLists [[0, 0, 4, 3, 2] | ||
| ,[0, 1, 3, 4, 3] | ||
| ,[1, 2, 3, 3, 4] | ||
| ,[0, 1, 2, 1, 0] | ||
| ,[0, 0, 1, 0, 0]] | ||
|
|
||
| -- | 'deconvolve' is used to correct the systematic error of blur (loss of | ||
| -- contrast in smaller features) in images. This particular variant is used | ||
| -- when no information about the distortion (blurring and noise) is known. | ||
| -- More info about deconvolution at <https://en.wikipedia.org/wiki/Deconvolution> | ||
| -- | ||
| -- <<images/frog.jpg>> <<images/frog_deconvolve.jpg>> | ||
| -- | ||
| -- Usage: | ||
| -- >>> u0 <- readImage' "images/frog.jpg" :: IO (Image VS RGB Double) | ||
| -- >>> let p = randomBlur | ||
| -- >>> let d = convolve Edge p u0 | ||
| -- >>> let us = iterate (deconvolve p d) d -- Iterative deconvolution | ||
| -- >>> u1 = us !! 1 -- one iteration | ||
| -- >>> u2 = us !! 20 -- twenty iterations | ||
| -- >>> let output = makeImage (rows u0, cols u0 * 4) | ||
|
Owner
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. This is an extremely inefficient way to concat images together, try using
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @lehins Sure sir, sorry I was unaware of that. Will change it.
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @lehins Sir, I tested the code using
Owner
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more.
Without some numbers/benchmarks I will not believe that. You probably aren't noticing the difference because deconvolution is an expensive operation and 20 iterations take long time. The operation of concatenation of 4 images the way you have it implemented in the example on its own will be much slower than if you use No worries on the delay. There isn't a particular rush on any of this, especially considering that GSoC is over. Congrats on a successful project. ;)
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @lehins I think I might be doing something wrong. I implemented it like :
Owner
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. You are doing it right.
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @lehins Very sorry, I actually did mean runtime but wrote compile time instead. Okay, I'll try to record the time. |
||
| -- >>> (\(r,c) -> | ||
| -- >>> let (i, c') = c `quotRem` cols u0 | ||
| -- >>> in index ([u0,d,u1,u2] !! i) (r,c')) | ||
| -- >>> :: Image VS RGB Double | ||
| -- >>> writeImage "images/frog_deconvolve.jpg" output | ||
| -- | ||
| deconvolve :: Image VS X Double | ||
| -> Image VS RGB Double | ||
| -> Image VS RGB Double | ||
| -> Image VS RGB Double | ||
| deconvolve p d u | ||
|
Owner
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. This algorithm isn't restricted to RGB images only. Also restricting to
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @lehins Okay sir, just one thing. I'll try my best to address the feedback on all my PR's asap but I've lectures and lab sessions today so won't be able to work on them before night time. Hope you understand. Thank you!
Author
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. @lehins Sorry sir, I didn't exactly understand. Do you mean that I should take input images as
Owner
There was a problem hiding this comment. Choose a reason for hiding this commentThe reason will be displayed to describe this comment to others. Learn more. No, what I meant is, it should be more polymorphic. So, try to fix up the code so the type signature is: |
||
| = u * conv (transpose p) (d / conv p u) | ||
| where conv = convolve Edge | ||
|
|
||
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Anything that is needed for testing, suppose to live in the test suite
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@lehins Okay sir, so this should be in the
usagesection in the documentation?There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yes you could move it to function documentation as an example. In that case I'd also recommend collapsing the big example (note
====__Examples__):In order to check out the affect, run
stack haddock --openThere was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@lehins Okay sir, will change it.