Gradient HDR compression, take 2

[georg.viehoever]

Hi,

this is a new version of the module I published in http://pixinsight.com/forum/index.php?topic=3114.0. As before, this module is available as source only, and has been tested on Fedora14-x64 only, but I don't expect any issues with compiling it on different platforms. It modifies an image by recomputing its gradients. For the purpose of HDR compression, we remove or reduce any large gradients. While Carlo's HDR module appears to be optimized for natural images, the approach used by this module is geared towards astromical images, preserving as much detail as possible (even if this results in an unnatural look). I will show some examples in later posts.

This version I would consider a beta version. Funtionality is largely complete and stable. Missing bits involve progress and cancel, and maybe some error handling.
Changes compared to the previous version:

- I switched from a sparse linear system solver to an FFT based solver. This is much faster than the previous version (19 seconds for a 10 MPixel image on a slow virtual Linux machine). It is also more accurate than popular alternatives such as multigrid methods.
- The current version uses the FFTW library (http://www.fftw.org/), which comes with Fedora 14, is available for all relevant platforms, and fast. Unfortunately, its license does not allow for distribution with commercial products. I will switch to a PI based version as soon as the necessary functions are available with similar speed within PI.
- This version is fast enough to allow for a realtime preview (!), so you can easily experiment with the sliders.
- The code has largely been ported to PCL based data structures. That avoids a number of copies, and the corresponding waste of memory. We also no longer need the GMM+ library.
- I have implemented a special mode that tries to preserve the colors of the stretched images ("Preserve Color" checkbox), by keeping the relative proportions of R:G:B unchanged.

I would love to hear if you find the module useful.

Georg

[georg.viehoever]

The attached image shows M42 (thanks to Harry) , the GradientHDR GUI, and the preview. In this example, we used the max gradient slider to limit the maximum gradients used. You can see how much details the process reveals in dark as well as in bright regions.

A problem of the process is that it can create a global gradient within the image. The reason: assume you have a 1D image with values 0.0, 0.5, 1.0, 0.0. The gradients are  0.5, 0.5 and -1.0. If you clip all gradients to 0.5, the gradiets are 0.5, 0.5, -0.5. The reconstructed image therefore has values of 0,0.5, 1.0, 0.5. You may need to use DBE or ABE to adjust for this.

Anyway, I don't expect the images generated by the process to be used directly. Most likely, they will be combined with other image, enhancing their level of detail.

Georg

[Harry page]

Hi

All looks very exciting , can't wait for a windows version  8)

Harry

[georg.viehoever]

Final example: A HDR photo of the window in our kitchen. I think this natural photo is quite useful for seeing the effect of this module. All pictures are treated with STF only.

- left: Original HDR shot, combination of 4 exposures, combined using the HDRComposition process. It is impossible to see details of the building in the background, the plant leaves and details of the tiles at the same time, no matter how you tune the STF.
- middle: GradientHDRCompression with "Preserve Color". Notice that you see details of the building in the background, the color variation in the leaves, and the structure of the tiles. In the area marked with "Preview 1" you can see a marker left by the builders of our kitchen, a "52". We never noticed, and we are living here since 6 years now. Even if you know it is difficult to find with the bare eyes.
- right: Same without "Preserve Color". Much more pale. It is possible to restore some of the color using the ColorSaturation tool, but never with the amount of detail visible in the other version.

Enjoy,
Georg

[Carlos Milovic]

Hi Georg!
The results looks very nice! Good job (could you upload those HDR images to test other algorithms? I'll send to you some links to other images that I use).

Did you see those gradients with the linear solver? Try a 3px zero padding of the image. That prevented artifacts in my DCT implementation.
BTW, I'll rewrite the solver today, using PI's FFT to calculate the DCT.

[georg.viehoever]

Hi Carlos,

thanks for the kind words. The kitchen image is available if you tell me where to upload it. Regarding M42, you need to ask Harry for the two images that I combined.

Yes, those gradients are visible with the linear solver as well. They are not related to the border, but are obviously caused by the imbalance of gradients that exists after truncating things.  I think your HDR implementation (based on http://www.cs.huji.ac.il/~danix/hdr/hdrc.pdf ?) avoids this problem by combining wavelets into the equation.

If you want to implement the PI based DCT into this module: its in gradients.cpp...

Georg

[Harry page]

Hi Carlos

The M42 HDR stack is on my FTP server
 
If you want it I will send you details on how to get it

Harry

[Carlos Milovic]

Hi Georg
Since we need to "fine tune" the processes to astro images, I think that I'll need just Harry's images.
For daylight purposes, I'm using mainly the radiance maps provided here:
http://www.cs.huji.ac.il/~danix/hdr/results.html
specially the Belgium house.

In fact, the wavelets approach was my idea They solved that using increasing radius gaussians and increasing size finite diferences. My calculation is a bit different, but yields similar results.

Right now I'm working on the FCT over my module... from nearly 1 minute on the belgium house (~1Mpx) it came down to less than a second. I'm having some artefacts and problems with aliasing... I have an idea how to fix that, but it will slow down the solver a bit. I'll let you know if I have sucess.
Ah, btw, I'm not using 2 dimensional DCTs... I'm solving the trilinear matrix that appears after 1D row DCTs. I think that this should be faster than another DCT pass (and another inverse DCT), since I'm just doing 2 O(n) operations...


Harry, I'll send to you a PM.

[georg.viehoever]

Ah, btw, I'm not using 2 dimensional DCTs... I'm solving the trilinear matrix that appears after 1D row DCTs. I think that this should be faster than another DCT pass (and another inverse DCT), since I'm just doing 2 O(n) operations...
..

Wow, O(n). Let me know when the solver is ready for showtime.
Georg

[Carlos Milovic]

I still have the O(w*log(w)) operations from the single DCT... here is the current solver. As you'll see, is pretty straightforward to include in any PCL based code, since it only needs the input image.

There is one problem I can't figure out how to solve: in front of very sharp edges the solver tends to generate a "gibbs effect" . Since I'm zero padding the borders in my HDRC implementation, this is generating those artefacts around the right and left borders... I'm looking for some workarounds, but doesn't seem to be easy

[georg.viehoever]

Hi Carlos,

can you give me an image that produces such a gibbs effect. Its not clear to me why a solver should show this effect.

Georg

[Carlos Milovic]

Forget about that I already fixed it, changing the algorithm... Now it takes nearly twice the time, since I had to double the length of the rows.... small price for VERY accurate results.

Here is the new solver.

[Carlos Milovic]

Another solver... using 2D DCT... to my surprise, elapsed times are almost the same, if not a bit faster in this case. Maybe it is due the optimizations the FFT transform has in the PCL, against my non parallel routines.


Georg, now I'm experimenting quite a bit with the HDRC algorithm... and I think that I'm seeing the same kind of large scale gradients you do, when I use strong gradient attenuations. I remember seeing that somewhere...

[Carlos Milovic]

Hey Georg

I uploaded a new source code of my module. Take a look at the DePadding function. There I'm performing a 2D linear fit to the gradient in the "border". I had previously done a zero padding there, so I was sure that it was a uniform zone. Then, subtracting the linear gradient did a wonderful job (I'm working in the log space, so it is a multiplication in normal space).

I'm going to read now your code, and see if I can modify it to use my solver... Once it is done, could you perform some execution time tests?

My module, with the whole process activated, took ~3secs to process the Belgium House (http://www.cs.huji.ac.il/~danix/hdr/pages/belgium.html).

[Carlos Milovic]

I think that all you have to do is include the solver header (solver_dct2.h or solver_dct3.h), and then change the line 322:
from: solveImage(li,rResultImage_p);
to:
__SolvePoisson( li );
rResultImage_p = li;