Night sky image processing - Part 7: Automatic Star Recognition

In the parts 1-6 of my "Night sky image processing" Series I wrote about the different stages of night sky image processing. In this article I want to put all the pieces together and develop an "automatic star recognizer" which takes an astro-image as input and outputs a list of the recognized stars with their HFD and FWHM values.

Basically the processing pipeline is as follows:
ROI = Region of interest

In order to achieve this I use most of the concepts I have looked into earlier: As input I used the following FITS image for testing. Continue reading →

Night sky image processing - Part 5: Measuring FWHM of a star using curve fitting - A Simple C++ implementation

In Part 4 of my "Night sky image processing" Series I wrote about star centroid determination with sub-pixel accuracy. Based on this centroid position the FWHM determination of the star takes place. The FWHM (Full Width Half Maximum) value is a measurement for the star width. In this part I write about the so called curve-fitting which is used to determine the FWHM value from such an image.
Note that for the following calculation and implementation the sub-pixel determination of the centroid is not considered. Let's say we want to determine the FWHM in x-direction (red line). There is of course also one in y-direction. Basically, what happens is:

  • The pixel lines through the centroid in x direction are extracted (those gray-level values usually form a Gaussian distribution - a "data-point" is defined as (x,y)=(pixel x-position, pixel gray-level value).
  • Based on those (x,y) data-points the 4 parameters of a Gaussian curve are determined so that the Gaussian curve and the data-points fit "as good as possible" - the algorithm used is the so called "Levenberg-Marquart" algorithm which tries to minimize the quadratic error between the data-points and the curve.
  • The result is a Gaussian curve i.e. a parameter set which describes the curve (c=center, p=peak, b=base and the w=mean width). One of those parameters - the mean width is the FWHM value.
The C++ implementation below uses the GSL (GNU scientific library) to make use of the Levenberg-Marquart algorithm. The primary application of the Levenberg–Marquardt algorithm is the least squares curve fitting problem. Imagine you have a set of data points Continue reading →

Night sky image processing - Part 4: Calculate the star centroid with sub-pixel accuracy

In Part 3 of my "Night sky image processing" Series I wrote about star-clustering. The result of this step was a list of stars and their pixels.

The next step is to find the center of each star. In this example the algorithm is only applied for one single star which is loaded from a FITS file for the sake of simplicity. The algorithm implemented here was proposed in "Improving night sky star image processing algorithm for star sensors" from Mohammad Vali Arbabmir et. al..

The implementation could be further optimized but to me the code below is a good compromise between clarity and efficiency. The paper mentioned above will probably help a lot to get a better understanding of the implementation. Generally there are two steps:

  • Step 1: Calculation of the intensity weighted center (IWC) (to get an idea of where star center might be) The calculation of IWC is based on the center of gravity CoG. The CoG is the same calculation as in physics but in this context just applied to an image. Each pixel brightness value has a weight - the brighter the "heavier". For example to get the "center" x position (this is the x position where the image has in x-direction the "brightest" value in mean. "In mean" means that all pixel values are considered). This directly leads to the following formulas for and :

    The only thing that changes between the two formulas is Continue reading →

Night sky image processing - Part 3: Star clustering

featured_image_fb_star_clustering_1200x630This article is about finding all the pixels which belong to a star - also known as star clustering. A short C++ implementation based on CImg is shown below. In Part 2 of my "Night sky image processing" Series I was writing about thresholding. The result of the thresholding is a binary image. In this context it means the "1" pixels belong to stars, the "0" pixels belong to the background. I used Otsu's method to calculate an according threshold value.

The binary image is now further examined in a process called clustering. In simple words the pixels which belong together (i.e. the neighbors) are grouped together. The result of this clustering process is a list of pixel groups which belong together - i.e. a list of stars and the pixels belonging to each star.

The implemented algorithm is based on the paper "Improving night sky star image processing algorithm for star sensors" from Mohammad Vali Arbabmir et. al. However, in detail it is slightly different. My implementation assumes that the supplied image cannot be modified during the process. Hence I am using a set to put all white pixels into. To me this is a good compromise between memory consumption, performance and simplicity. However, the implementation can be further improved with respect to runtime performance when one would used a 2D array for the white pixels instead of a set.

The clustering algorithm is relatively straight forward:
  1. Collect all white pixels of the binary image to a bin
  2. Continue reading →

Night sky image processing - Part 2: Image binarization using the Otsu thresholding algorithm

This article is about image binarization using the Otsu method. A simple C++ implementation based on the CImg library is shown. In Part 1 of my "Night sky image processing" series I have taken a look on anisotropic diffusion in order to reduce the image noise.

The next topic I want to have a look at is thresholding. In this step it is decided which pixels belong to the background and which pixels belong to a star. Before it can be decided which pixels belong to a star (also called clustering - see part 3) it first has to be generally decided which pixels belongs to any star and which pixels belong to the background. This step is called thresholding i.e. converting a grey-scale image into a binary image. There are many different thresholding algorithms out there. All have there advantages and drawbacks. I found ImageJ very helpful to test different thresholding algorithms. After trying different ones I found that Otsu's method does a quite good job. After further reading I found out that this method is especially useful for images with bi-modal histograms. To my knowledge night shy images usually do not have bi-modal histograms. However, so far Otsu's method still gave very good results with my star images. Furthermore it is relatively efficient and simple to implement. For those reasons I decided to go this way.

Otsu's algorithm tries to find the threshold so that the variance of all pixel values in the foreground and the variance of all pixel values in the background becomes minimal. The variance is a measure of the spread. The smaller the variance of for example the foreground pixel values, the closer the gray levels of the pixels and the more likely it is that they belong together. At least that's how I understand Otsu's method.

The following formula calculates a kind of "mean" variance of the foreground variance and the background variance . In fact it is no mean but a weighted sum of both. The weights are the number of pixels in the respective group:

while is Continue reading →