Skip to content

FEAT: Image Deconvolution Functions #1881

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

Merged
merged 1 commit into from
Jul 13, 2018
Merged

Conversation

9prady9
Copy link
Member

@9prady9 9prady9 commented Jul 21, 2017

Iterative deconvolution algorithms

  • Landweber
  • Richardson-Lucy

Linear deconvolution algorithms
- [ ] Wiener

  • Tikhonov

Wiener has been disabled for now. I haven't had a break finding the bug in that. We can enable the API in future.

@9prady9 9prady9 added this to the v3.6.0 milestone Jul 21, 2017
auto ONES = createValueArray(dims, scalar<T>(1));
auto HALF = createValueArray(dims, scalar<T>(0.5));
auto midVal = arithOp<T, af_mul_t>(mxArr, HALF, dims); // B/2
auto diff = arithOp<T, af_sub_t>(prev, midVal, dims); // |I(n-1) - B/2|
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You are not calling abs anywhere.

Copy link
Member Author

@9prady9 9prady9 Jul 21, 2017

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

missed that. Need to look into why tests passed, or if the abs is needed or not. Thanks. I know why they passed - this(jansson) algorithm is not being tested right now as I don't have reference data. ITK has Landweber and I added tests only for that as of now.

@9prady9 9prady9 changed the title [WIP] FEAT: Deconvolution Functions [WIP] FEAT: Image Deconvolution Functions Jul 21, 2017
@9prady9 9prady9 force-pushed the deconvolutions branch 5 times, most recently from 1c29b6c to a0274c1 Compare July 27, 2017 09:35
@pavanky pavanky changed the base branch from devel to master August 23, 2017 15:21
@9prady9 9prady9 force-pushed the deconvolutions branch 2 times, most recently from 4cc2163 to 5c25132 Compare September 12, 2017 11:40
@9prady9 9prady9 force-pushed the deconvolutions branch 2 times, most recently from d1eec82 to 94a167c Compare December 13, 2017 03:37
@9prady9 9prady9 changed the title [WIP] FEAT: Image Deconvolution Functions FEAT: Image Deconvolution Functions Dec 28, 2017
@9prady9 9prady9 force-pushed the deconvolutions branch 2 times, most recently from 815b1cc to d851935 Compare January 25, 2018 05:10
@9prady9 9prady9 force-pushed the deconvolutions branch 2 times, most recently from 468351c to 0f04d60 Compare February 13, 2018 08:18
@9prady9 9prady9 requested review from umar456 and pavanky February 13, 2018 09:53
@9prady9 9prady9 dismissed umar456’s stale review February 22, 2018 14:29

New changes pushed.

@9prady9 9prady9 requested a review from umar456 February 22, 2018 14:29
@umar456 umar456 requested a review from syurkevi February 23, 2018 00:54
@umar456 umar456 removed this from the v3.6.0 milestone Feb 27, 2018
@9prady9 9prady9 removed the request for review from syurkevi April 3, 2018 04:03
umar456
umar456 previously requested changes Apr 8, 2018

Iterative Deconvolution Algorithms

The following table shows the iteration update equations of the respective deconvolution algorithms.
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Can you make this more descriptive? It doesn't talk about what the algorithm does.

Copy link
Member Author

@9prady9 9prady9 Apr 9, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Sure. I think few pictures can help illustrate the case too.

- \f$ \alpha \f$ is the relaxation factor
- \f$ \otimes \f$ indicates the convolution operator

\note The type of output \ref af::array from deconvolution will be double if the input array type is
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Format lines to be 80 characters per line?

@@ -13,6 +13,7 @@ project(ArrayFire-Example-Image-Processing

find_package(ArrayFire)

set(ASSETS_DIR "${ArrayFire_SOURCE_DIR}/assets")
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

You probably need to update to reflect our most recent changes to the example folder.

C++ Interface for Iterative deconvolution algorithm

\param[in] in is the blurred input image
\param[in] ker is the kernel(point spread function) known to have caused the blur in the system
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

80 characters per line.

/**
C++ Interface for Iterative deconvolution algorithm

\param[in] in is the blurred input image
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What type of images do you expect? Are there any restrictions?

Copy link
Member Author

@9prady9 9prady9 Apr 18, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Addressed. Check out the docs for this function to see the updates.


\ingroup image_func_iterative_deconv
*/
AFAPI array iterativeDeconv(const array& in, const array& ker,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I feel like the name should be deconvIterative and deconvInverse. It emphasizes the operation you are performing.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Having deconv as prefix sounds odd and reverse. Having Deconv as a suffix has an ease to the name when we want to talk about the function or refer to it somewhere. Also, it has a nature flow as to how the name translates to the brief description of the algorithm.

iterativeDeconv - Iterative Deconvolution
inverseDeconv - Inverse Deconvolution

sounds better compared to

deconvIterative - Iterative Deconvolution
deconvInverse - Inverse Deconvolution

// cuFFT/clFFT works with any data size, but is optimized
// for size decomposition with prime factors up to
// 7.
const dim_t GREATEST_PRIME_FACTOR = 7;
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Does this mean that the OpenCL and the CPU versions will behave differently? We should avoid doing that if that is the case.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

There is no difference in behavior except that I am using prime factor that each backend is optimized for when we use ffts. The error metric computed against the test data is same across all backends. The prime factor is used to determine the padding length, nothing else.

std::vector<OutType> goldData(nElems);
ASSERT_EQ(AF_SUCCESS, af_get_data_ptr((void*)goldData.data(), goldArray));

ASSERT_EQ(true, compareArraysRMSD(nElems, goldData.data(), outData.data(), 0.03));
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How did you get the 0.03 threshold for this operation? Wouldn't this be depended on the type of the input and output data and its range?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

I guess what I am asking is that have you done any comparison on the image itself and looked at the results at the pixel level and calculated why there are differences. Sometimes there are subtle bugs which cannot be viewed or appear when you look at metrics like RMSD

Copy link
Member Author

@9prady9 9prady9 Apr 18, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

It would depend on the input range unless we bring them into the same range, which I am doing both while generating the test data from ITK and after calculating the output in ArrayFire tests.

I have checked the images in GIMP and I didn't notice any odd/huge artifacts. If there is any difference, it is probably some(random spread out) pixels here and there I think.

namespace kernel
{
static dim_t
calcIdx4Border(const dim_t i, const dim_t lb,
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This function should really be optimized. This is being used a lot in your function and it is using some expensive operations (lots of branching and modulus). You should think more carefully about how to implement this.

- \f$ \gamma \f$ is a user defined noise variance constant


\note The type of output \ref af::array from deconvolution will be double if the input array type is
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What is the range of the float values. Is it 0-1 or 0-255?

Copy link
Member Author

@9prady9 9prady9 Apr 9, 2018

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

The algorithm/function itself is niether dependent on the range of the data nor alters the range of the data.

@9prady9 9prady9 force-pushed the deconvolutions branch 2 times, most recently from ab90678 to 64cfac8 Compare April 19, 2018 11:35
@9prady9 9prady9 requested a review from umar456 April 19, 2018 11:36
@9prady9 9prady9 dismissed umar456’s stale review April 19, 2018 11:37

Addressed feedback

@9prady9
Copy link
Member Author

9prady9 commented Jun 6, 2018

@umar456 I have cleaned the style a bit further. It is ready for merge unless there is more feedback.

The following iterative algorithms have been added
* Landweber - Tests PASS
* RichardsonLucy - Tests PASS

The following linear inverse algorithms have been added
* Tikhonov - Tests PASS

Wiener Filter  - Tests FAIL - Disabled temporarily, The function
is also not exposed at the moment.

An example demonstrating both iterative and inverse deconvolution
has been added as well.

Unit tests using grayscale images have also been added.
@@ -273,3 +273,5 @@ make_test(SRC where.cpp)
make_test(SRC wrap.cpp)
make_test(SRC write.cpp)
make_test(SRC ycbcr_rgb.cpp)
make_test(SRC inverse_deconv.cpp)
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Combine the tests into one file.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Will in another PR while doing cleanup with other functions.

richardsonLucy(currentEstimate, paddedIn, P, Pc,
iters, normFactor, odims);
break;
default:
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

This should throw an error if the wrong algorithm value is sent.

@9prady9 9prady9 merged commit d88f626 into arrayfire:master Jul 13, 2018
@9prady9 9prady9 deleted the deconvolutions branch July 13, 2018 07:00
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Development

Successfully merging this pull request may close these issues.

3 participants