Better "debayering" strategies for OSC cameras

[Ignacio]

I have been thinking about this for a while now, motivated to improve color resolution and color noise in and around small stars (that take a few pixels), although the arguments here may apply in general.

Is there a better way to take advantage of the spatial color information embedded in a series of raw bayered images that are slightly missaligned? I have observed that debayering low S/N subexposures can introduce quite a bit of artifacts and noise, even if done after calibration and with the better methods. Besides, debayering algorithms have been developed for higher S/N images and have a hard time dealing with steep, contrasty edges, such as we find in stars. Also, they don't take into account info from other subs in a session.

What if "debayering" is done at the time of stacking, as we do by mixing RGB separate images taken with monochrome cameras? For the sake of argument, take the case of a star that illuminates a single pixel. In one sub, is a red pixel, in the next sub it illuminates a green pixel, and so on and so forth. For the complete (dithered) session, for this particular star, we will very likely have sample statistics for the three colors. Think of it as randomly turning an RGB filter wheel, but at a pixel level. Granted that for a given spatial position (ie, star), you may end up with more data for one color channel than the other. But it seems to me that the color info in the stack is "pure" raw data, that has not been altered by debayering interpolation. One would end up with, say, an R channel image, where each pixel has a different amount of integration time, which reflects physical reality. From there one can do all sorts of processing (alla interpolation), but with a much better signal-to-noise starting point.

Does this make sense? Can this be tried somehow with the current tools in PI 1.7? I would like to experiment.

Thanks,
Ignacio


     

[Carlos Milovic]

Yes. I believe Super-resolution is the approach we have to follow here.
Roughly, it states that there is an ideal image, of twice (or more) the resolution at each axis, and all we have are subsampled, misaligned and blured versions of it. By solving the inverse problem, with a given regularization, it is possible to find that ideal image. I've got a paper on Wavelet-based superresolution that was used for astronomical images that may be well suited in this case. Also, I'm working with an state of the art regulatization, that may be incorporated in the scheme. All of this is in our development plans for the next year.

[Ignacio]

Sounds interesting, Carlos. I would like to read more about this. Can you point me in the right direction?

Also, is there a simple version of this that we could apply with current tolls?

Thanks again,
Ignacio

[Carlos Milovic]

This is the paper I was talking about:
www.ee.duke.edu/~willett/papers/WillettADASS03.pdf

I already implemented the registration algorithm that is suggested there, but it still needs work (lacks some robustness and it hasn't being optimized enough; for small images is pretty fast, but for large images StarAlign is faster and also more accurate). Implementing the EM steps with the starlet transform should be quite straightforward with the PCL. Using the new regularization scheme will require more development.

[Ignacio]

Thanks, Carlos, I will take a look.

In any case, how does image registration work when applied to a set of CFA RGB images? That is, what happens with the "dark" pixels from any of the color channels when they get shifted a non-integer number of pixels? Are they ignored or do they carry a zero value to the new array location?

Ignacio

[Carlos Milovic]

Well, I think that for registration purposes, they should be debayered first (any method), to get he most accurate displacements and rotations of the whole set. Integration is done with bayered data. So, as you may suspect, this is not a fast procedure. Even less if you consider a deconvolution involved.

To be more accurate, what happens is that we calculate an ideal image (we may use a standard integration from all debayered subframes as a first  approximation, 2x upsampled) and then it is compared with the actual undersampled frames that we captured by traslation/rotation, bluring and finally undersampling and multiplying with the bayer matrix pattern. So, we end with a data fidelity term that compares the projections of the ideal image to each of the images we captured. Then comes the regularization step.... And everything repeats until convergence.

[Ignacio]

On your first point, I know that standard registration requires debayered frames. Still, I would like to understand how the data is treated. I would guess that the translation/rotation matrix is applied to the pixel array of each color channel, right? If this is so, pixels that due to the bayer matrix have zero value, would carry thru the transformation and (partially) land in another pixel, averaging down it's value. Right? If so, it would preserve the integrity of the data on each channel, regardless of the accuracy of the registration.

In your second paragraph you lost me, Carlos.   I'll have to read the paper you mentioned first, I guess.

Thank you,
Ignacio

[Carlos Milovic]

No, the process is done the other way around. Instead of registering all the subframes to the reference frame, you register the reference frame to all the subframes.
Let's elaborate a bit more. We start with the set of images (unregistered, CFA). We select one of them as the reference, and then the register all the images, after debayerization, as usual, but instead of applying the transformation, we just save the transformation matrices (or, we may transform them, and average to get a new reference image, that has higher SNR). Now, we create a first approximation to the ideal image. It could just be the reference image, upsampled 2 times. The goal on the inverse problem is to minimize certain potential, that contents two terms: a data fidelity term, and a regularization.
The data fidelity term is the total absolute diference between the ideal image to each of the captured frames. So, to compute that term, we have to translate/rotate the ideal image (RGB, normal color), then undersample it (integer resample, averaging), and finally subsample it according to the CFA pattern. We end with several greyscale, bayered images that are the projections of the ideal image to each captured RAW CFA frame. Then we calculate the difference.
The regularization step may be understood just as a noise reduction filter. In fact, the expectation-maximization algorithm is just these two steps (data fidelity and then noise reduction). You alternate both steps until a convergence is reached, and the optimal image is achieved, given our potential terms.

[Ignacio]

Ah! got it, very interesting. In this way you do not disrupt the bayer array of each subexposure, but you need to build a reference bayered frame for each one, to compute some measure of distance. Clever.

I look forward to your implementation, keep us posted.

best
Ignacio