Help required with BM3D Image denoising algorithm implementation

Question

This is what I am doing and will be grateful for any help.

1. I selected 16 similar blocks to a reference block from a search block, put them in a 3D array, Wavelet Thresholded them, took inverse transform.

2. These 16 similar blocks belong to one 32 * 32 search block as mentioned above.

3. Now I should multiply each of these blocks with their weight, and then add them all up, creating a 2D array from a 3D array.

4. This new block then replaces all the 16 original blocks that I had selected before for hard thresholding.

Am I correct here? I will be absolutely grateful for a response.

Answer 1

Step 3 is not correct. The first part is ok (each denoised result is weighted by some likelihood of correctness, that was obtained previously from the variance of the data), but you don't have to add them up. Each slice from the 3D array is a 2D patch that is back-projected to the resulting image.

The sum happens when 2 back-projected patches do overlap in the final image domain.

About back-projection

Back-projection is the process of taking a 2D slice from the denoised stack and put it back in image space. This is required because unlike most variants of NL-means BM3D denoises patches and not the central pixel of a patch.

So, you have 2 images when implementing BM3D:

• the original (noisy) one;
• the result (denoised) image.

The spatial location of a patch should not change during the denoising process, otherwise the image content would be destroyed. So, each slice from the denoised patch stack goes back to its original place.

Now, one issue remains: what to do when back-projected patches overlap? What value should we keep for a single pixel? The solution chosen by BM3D is that each pixel contains a weighted mean of the patches that overlap in each pix.

To implement this, you can simply pre-multiply each slice of the denoised stack by its weight, add a third image to store the sum of the weights inside each pixel, sum for each pixel its value from all the denoised patches that contain it, and in a final step (once all the patch stacks haves been processed) divide the result image by the weight image.

Image processing main loop (outline)

There is no "unleft" pixel in the final result. The outline of the main loop is as follows:

1. for every pixel of the image, form a reference patch;
2. for a given reference patch, find n <= N (B=16 in your case) patches whose distance to the reference patch is less than some threshold and make a stack;
3. denoise the stack;
4. back project the patches in the denoised stack;
5. back to 1) for the next pixel).

Thus, you always have one stack for a given patch (degenerate case: only 1 patch in the stack), and a pixel is always involved in several patches.