Here's an example image to make the problem very clear:
I've used OpenCV's blob_log() as well as the rank.otsu threshold approaches. The problem with the blob detection is that the blob radii are quantized, e.g., 3, 7, 11, etc. The otsu local adaptive threshold works fairly well but needs severe hyperparameter tuning -- in other words, it works in some locations but not others.
I don't know how I can ask the following question: "find all circular blobs in the image. If the blob is large (meaning it is overexposed), expand the search to include diffraction spikes (at known fixed angles) plus the Airy disk."
My goal is to remove the stars and their optical effects from the image, then use an inpainting algorithm to fill in the missing pixels. I have several inpainting candidates to try.
There is noise in the image, and I'm fine with missing any stars that are indistinguishable from noise.
This is one channel -- the other image channels do add additional information, and the actual removal of stars/effects will be based on the logical or of all three analyzed channels.
I am aware that there are other methods for star removal that pop up in google search. The best automatic one I've used is starnet, e.g. in https://github.com/nekitmm/starnet . The issue I have with it is it doesn't separate out the star mask but goes ahead and does the inpainting -- unfortunately, very poorly in large regions. And it doesn't handle spikes well, per the author's commentary and my observations.
I wanted to avoid a machine learning approach. Instead, I wanted to solve the problem directly from first principles. It seems feasible to build a structural model of the star/spike/airy disk with some effort(e.g., as central disk at max intensity becomes larger, spikes grow larger and so does Airy disk.) But how do I match this model to the objects in the image?
Thanks in advance.
The desired output is a 1 bit mask image. Bit set to 1 if this is a star pixel, 0 otherwise. I don't have a problem converting a list of (y,x,r) values to a bit mask, and already have code to do that.
The blob detectors I have used (a) don't compute accurate blob sizes (b) generate blobs for some of the noisy areas.
Regarding finding the isolated blobs being a snap, give me an idea of what you consider is a "snap." My very first attempt looked for local maxima and attempted to grow a circular region. It was very slow and had problems with noise.
The region-growing approaches all failed miserably. A different 2d filtering approach turned out to be more promising. The basic idea is to:
- Capture all of the overly bright areas via thresholding, e.g., against around 248 in the bottom image. That gives the first mask = image > 248, e.g.
- Apply an emboss filter at 45 and 135 degree angles (filter kernels are [-1 -1 0][-1 0 1][0 1 1] plus 128 for the 45 degree angle.)
- Threshold values more than 45 away from grey (128). mask = (embossed < 128-45) | (embossed > 128 + 45), e.g.
- Or together all of the three thresholds.
That gives the following result:
That approach has two magic numbers, yes. But it works fairly well for such a simple algorithm.
I'm aware the blobs are a bit too big in places. Also, the blue and red channels of the original image don't work quite as well due to more noise in them. Time to experiment further.