A&A 457, 729-736 (2006)
DOI: 10.1051/0004-6361:20064852

Optimal filtering of solar images using soft morphological processing techniques

S. Marshall¹ - L. Fletcher² - K. Hough¹

1 - Department of Electronic and Electrical Engineering, University of Strathclyde, Royal College Building, 204 George Street, Glasgow G1 1XW, Scotland, UK
2 - Department of Physics and Astronomy, University of Glasgow, Glasgow G12 8QQ, Scotland, UK

Received 13 January 2006 / Accepted 27 June 2006

Abstract
Context. CCD images obtained by space-based astronomy and solar physics are frequently spoiled by galactic and solar cosmic rays, and particles in the Earth's radiation belt, which produces an overlaid, often saturated, speckle.
Aims. We describe the development and application of a new image-processing technique for the removal of this noise source, and apply it to SOHO/LASCO coronagraph images.
Methods. We employ soft morphological filters, a branch of non-linear image processing originating from the field of mathematical morphology, which are particularly effective for noise removal.
Results. The soft morphological filters result in a significant improvement in image quality, and perform significantly better than other currently existing methods based on frame comparison, thresholding, or simple morphologies.
Conclusions. This is a promising and adaptable technique that should be extendable to other space-based solar and astronomy datasets.

Key words: methods: data analysis - techniques: image processing - Sun: corona - Sun: coronal mass ejections (CMEs)

1 Introduction

The charge couple device (CCD) is the workhorse detector of many space-based astronomy and solar physics imaging instruments, such as the Hubble Space Telescope Wide-Field/Planetary Camera, and the Solar and Heliospheric Observatory (SOHO, Domingo et al. 1995) Extreme Ultraviolet Imaging Telescope (Delaboudinière et al. 1995) and Large Angle Solar Coronagraph (LASCO, Brueckner et al. 1995). Its use is also planned on many upcoming missions. CCD detectors are fairly robust, can be adapted for use over a wide range of wavelengths from X-ray to infra-red, and produce digital output that can be readily telemetered to Earth. These factors have made them ideal for use in space. The first orbital solar observatory using CCD technology (in its Soft X-ray Telescope) was the Yohkoh satellite (Ogawara et al. 1991) and since that time, CCD-based observational study of the solar atmosphere has proceeded in leaps and bounds.

In the space environment, CCDs are bombarded by cosmic ray ions (CRs) which register as counts in the CCD. Therefore for space-based use, CCDs are radiation hardened (which minimises permanent damage) and shielded. However, shielding leads to increased costs, and can never prevent energetic particles entering the telescope aperture. In particularly high-radiation environments, radiation hits can still significantly degrade the data. Such environments include; outside the magnetosphere of the Earth, particularly during active space weather conditions; low-altitude and high-inclination orbits; during passage of the South Atlantic Anomaly. CR hits can be problematic even for ground-based observations with long integration times.

In this paper we focus on imaging data taken with LASCO, onboard the SOHO satellite. This suffers significant particle hits (compare the left and right hand panels of Fig. 1), particularly during the early phases of coronal mass ejections (CMEs). High energy particles accelerated close to the Sun (Kahler 1992) - possibly at the flare acceleration site itself, or in the early stages of a CME-driven shock - impact the CCD leading to the characteristic "snow'' which in some cases completely obscures the very earliest phases of the dynamic event itself. Similar problems are seen in data from the Transition Region and Coronal Explorer satellite (TRACE, Handy et al. 1999) which is in a low-altitude, high-inclination orbit, taking it through the auroral zone of the Earth's radiation belts. It too suffers from particle hits, which are particularly troublesome in images taken in the extreme ultraviolet where counts are low.

	Figure 1: Two images from the LASCO C2 telescope, showing the effect of cosmic ray hits. The left hand image is relatively unaffected, whereas in the right hand image much of the interesting structure is obscured by "snow''.
Open with DEXTER

Various post-processing techniques are already available, and are used with success for CR noise removal in cases where CR noise does not too greatly affect the images. However they are of limited used in data badly corrupted by CRs. The present paper demonstrates the application of non-linear filtering techniques to astronomical and solar images, specially adapted for the removal of cosmic ray noise, and successful also in cleaning very noisy data. In the next section we give a brief introduction to some of the general approaches already used for CR noise correction, as well as the particular case for the LASCO images which we study here. This is followed in Sect. 3 by a description of morphological filters, the design method in Sect. 4 and the results in Sect. 5.

  
2 Cosmic ray noise correction

The spoiling of images by CR hits reduces the effective duty cycle of an instrument, and hence the return on investment. In the case of solar activity observations it can also obscure the crucial early phase of the event. Post-processing of images is often used to attempt to remove the effects of unwanted CR hits. The best post-processing method depends on the type of the observation. In sequences of images with a relatively low CR hit rate, and in which the target objects do not change significantly then CR pixels can be relatively easily flagged by co-aligning and comparing subsequent images (e.g. Shaw & Horne 1992). However, this type of technique often does not perform well in data where the image changes from frame to frame, as is the case in solar observations. In such cases, an algorithm capable of identifying and removing/replacing CR artefacts from single images, based on the known or assumed properties of these features, is likely to perform better. A frequently-used method is the median filter, in which pixels having counts above some defined threshold are replaced by the median value of some chosen number of neighbouring pixels. Three more sophisticated algorithms, tuned for the removal of CR spikes from stellar fields, are reviewed in Farage & Pimbblet (2005). These were based on identifying properties of CR artefacts: (i) their "sharpness'' (i.e. narrower than the instrument point-spread function); (ii) their non-Gaussianity (i.e. due to the instrumental PSF the sub-area around a stellar feature will have a roughly Gaussian distribution of counts per pixel above "empty sky'' background, from which high-count CR pixels can be separated); and (iii) their intensity above some user-defined threshold in spatially smoothed images. However, in the solar case, these algorithms are of limited use. For example we know that solar flares can produce saturated pixels in the UV and EUV on the scale of the instrument point-spread function, which change from image to image, and which are scientifically informative. Also, in a complicated solar field there is often no easy way to define a "background'' which does not risk losing physically useful information. Simple thresholding algorithms, often coupled with a comparison of saturated pixel locations in previous and following images, are in fact used for CR detection in solar cases. But it is certainly not an optimum method, failing to distinguish "real'' saturated pixels (i.e. due to processes on the Sun) from artefacts. Frame comparison techniques are also used with fairly good results (we show an example in this paper) but fail when the CR flux becomes extreme.

"Training'' methods, in which an algorithm (e.g. a genetic algorithm) is taught to recognise the characteristics of a CR hit or any other feature may provide an alternative in the solar case. Murtagh (1992) describes one such case, in which several properties of a desired class of object are first specified, producing a "filter'', and objects identified based on the (degree of) occurrence of these properties in a training phase and then assessed in a test phase. In this paper we describe and show the results of a case study in which a type of filter known as a soft morphological filter is developed and applied to data from the LASCO telescope onboard SOHO.

LASCO is a set of three nested coronagraphs, imaging the solar corona from 1.5 to 30 solar radii. The outer two coronagraphs are broad-band white-light instruments, imaged onto $1024 \times 1024$ pixel CCDs. The figures in this article show images from the C2 coronagraph, the field of view of which extends from 2 to 6 solar radii. Before transmission to ground, the data are transformed and/or compressed to reduce the telemetry requirement. The use of a lossy compression algorithm will also introduce spurious features into the uncompressed images, although at a much lower relative intensity than CR hits. As mentioned above, LASCO images can suffer badly from CR hits during a solar energetic event such as a coronal mass ejection, thus obscuring the very dynamics one hopes to study.

Figure 2: Upper: graphical representations of standard morphological erosion (LHS) and dilation (RHS) operations on a 1D signal by a circular stucturing element. The solid line is the signal received, and the signal following the operation is given by the dashed line. Lower: graphical representations of soft morphological erosion (LHS) and dilation (RHS) operations on a 1D signal by a circular stucturing element with a "hard'' centre and a "soft'' surround. The solid line is the signal received, and the signal following the operation is given by the dashed line.

Open with DEXTER

  
3 Soft morphological filters

Soft morphological filters are a non-linear image processing method particularly effective for noise removal. They originate from the field of mathematical morphology which was first developed by Matheron (1975) and Serra (1982) to model the propagation of gases through porous rocks. Initially the work applied only to binary images, but it has since been extended to greyscale (Serra 1988), colour (Angulo & Serra 2003) and complete lattices (Heijmans 1994). The soft morphological filter was first introduced by Koskinen & Astola (1994). Design techniques for these filters based on genetic algorithms have been developed and can be found in Kraft et al. (1997), Harvey & Marshall (2000), Marshall et al. (2000) and Hamid et al. (2003). Existing applications include spatio- temporal filters for film archive restoration (Hamid et al. 2003). It is not possible to introduce here the entire field of mathematical morphology, so instead a brief overview will be given with specific details for the implementation of soft morphological filters.

Mathematical morphology involves the study of signals and images based on their shape. The classical signal and image processing techniques seek to simplify the input data by mapping it onto orthogonal spaces. This results in superposition theorems in which the data is decomposed into independent elements which can be processed separately without cross interference. Examples are Fourier transforms where data is separated into independent sinusoids, and principal components analysis (e.g. Anderson 2003) in which the data are reduced to orthogonal basis functions describing the dataset variances.

On the other hand morphological image processing uses a sub-image (known as a structuring element) which acts as a probe. The morphological operation known as an erosion results from sliding the probe under the data and testing if it "fits beneath'' the surface of the image. The processed image is produced by tracking the locus of a reference point on the structuring element as it passes over the data. The size and shape of the structuring element determines the nature of the output image. A corresponding operation, known as a dilation results from sliding the structuring element over the surface of the image. Multiple sequential operations may be used. For example, an erosion followed by dilation is known as an opening, a dilation followed by an erosion is known as a closing. An example of a standard erosion applied to a one-dimensional signal is given in Fig. 2 (upper left panel). The structuring element is circular with its reference point in the centre. It can be seen that the overall level of the signal has been reduced. Also peaks which are too narrow to accommodate the structuring element are removed. This is effective for removal of noise. Other areas of the processed signal follow the same general shape. In image processing applications the structuring element would be three dimensional, i.e a ball rather than a circle. An equivalent example of a standard dilation with the same structuring element is shown in Fig. 2 (upper right panel).

Soft morphology is slightly more subtle than standard morphology described above. It represents a hybrid method which combines the structuring element approach of standard morphology with the fine tuning afforded by rank selection. In soft morphology the structuring element consists of a hard centre, $\alpha$ and a soft surround, $\beta$. The hard centre behaves in a similar way to the structuring element in standard morphology. That is, for soft erosion the whole of the hard centre must be "beneath'' the signal surface. On the other hand only a proportion of the soft surround, controlled by a rank parameter, r must be beneath the signal surface. An example of soft morphological erosion and dilation are given in Fig. 2 (lower left and right panels respectively). A soft morphological erosion, indicated by  $\circleddash$, is therefore defined by three parameters, $\alpha$, $\beta$ and a rank parameter r. Its output when applied to an image or signal f is given by

$$f \; \circleddash \;\left[ \alpha,\beta, r \right](x) r_{\rm min}^{\rm th}\left\{r\;\lozenge\left( f\left(y+x\right)\right.\right.$$ =

$$-\left.\left.\alpha \left(y \right) \right)\cup \left(f\left(z+x\right)-\beta\left(z\right)\right) \right\}$$

for $y \in \alpha$ and $z \in \beta$, where $r_{\rm min}^{\rm th}$ denotes the $r^{\rm th}$ smallest element of the set in curly brackets, and $\lozenge$is a repetition operator, such that $r \lozenge b(y)$ indicates that b(y) is repeated r times. Similarly soft morphological dilation ( $\circledplus$) is defined as

$$f \; \circledplus \;\left[ \alpha,\beta, r \right](x) r_{\rm max}^{\rm th}\left\{r\;\lozenge\left( f\left(y-x\right)\nonumber\right.\right.$$ =

$$\left.\left.+\alpha \left(y \right) \right)\cup \left(f\left(z-x\right)+\beta\left(z\right)\right) \right\}.$$ (1)

These equations imply that there are two basic steps to the implementation of soft morphological filters. The first is adjustment of the input pixel values: in the case of soft erosion (dilation), the structuring element values are subtracted from (added to) the pixel values within the filter window. The second stage involves the selection of the $r{\rm th}$ ranked value: that is, the resulting values from the previous step are first rank ordered and then the $r{\rm th}$(largest or smallest) value is identified. The fact that the rank parameter is the same as the repetition value associated with the hard centre ensures that the hard centre always remains entirely below (erosion) or above (dilation) the image surface. By careful design of the structuring elements, soft morphological filters can be used to remove noise from images. In this way the filter models the inverse process and produces an optimum mapping from noisy to restored image. We note additionally that since erosion and dilation are dual operations (meaning that they perform equivalent operations on the positive and negative features of the image) they are equally capable of filtering both positive and negative spikes. For example eroding the foreground of an image is directly equivalent to dilating the background.

The soft morphological filter used in this work was designed using a genetic algorithm (e.g. Holland 1995). The non-linearities in the filter make it difficult to produce a deterministically-designed optimum solution. Instead an iterative search approach is used to develop and train an optimum filter, as described in the next section.

   
4 Filter design and training

4.1 Creation of a training set

The first step in removing the noise from any image is to understand the nature of the disturbance. In this case, the distortion is caused by cosmic rays hitting the CCD in the LASCO telescope and causing the cells to overload, producing an image that suffers from extreme "white out''. These CR hits need not result in completely saturated pixels, therefore a means of correction based on the overall appearance rather than just a pixel intensity is sought. A key point in soft morphological filters is that as well as containing the noise rejection properties of rank order filters they also combine information about the structural content (shape) of the images. It is vital to ensure that the noise model used in training is appropriate to that affecting the real images. For this example, clean images from the SOHO telescope were taken from the LASCO data archive at sohowww.nascom.nasa.gov. These images were cropped from $1024 \times 1024$ pixels to $150 \times 150$ pixels, to reduce the training time (although too small a training set will result in the filter being undertrained). White marks, similar in shape to those seen in Fig. 1 were added to create the training data. This was done manually through cut and paste. Figure 3 shows examples of training images used. With experience we find that using 10 cropped images produces well-trained filters.

	Figure 3: Examples of two pairs of training sets generated for creation of a soft morphological filter for cosmic ray removal. The left-hand panel shows "ideal'' clean LASCO images, and the right hand-panels have had cosmic ray "noise'' added by hand. In total, ten of these sets were used in the training process.
Open with DEXTER

4.2 Training

The next step is to prepare and execute the actual training process using a genetic algorithm, in which the filter parameters and constraints are represented by a binary "chromosome''. Thirty filter chromosomes are initially created using a pseudo-random number generator. The filters are applied to the test data and the fitness assessed according to two quality measures, described below. Thereafter, the best 15 filter chromosomes are kept and are subject to the GA techniques of crossover and mutation (Holland 1995) to create new chromosomes from which to generate new filters. This propagates desirable features of the successful filters, and allows the introduction of "new'' and potentially useful information to the filter chromosomes.

The filter constraints on this particular training run were as follows:

• the filter used was to be a $5 \times 5$ soft morphological filter;
• the weighting parameters were to be horizontally and vertically symmetrical so only one quarter of the mask was being trained;
• there would be a total of 4 possible morphological operations per iteration with the possible operations being either erosions, dilations or no operations. For example, it would be possible for one of the sets of operations to consist of an erosion, a dilation, a no operation followed by another erosion;
• there were to be 30 chromosomes per iteration;
• the mutation factor was to be set at 0.03 probability of occurring and the crossover rate was 0.7 probability of occurring.

The fitness is measured by comparing the noisy and clean images with a quality measure and then the filtered and clean images with the same quality measure. Two training runs were performed using two different measures of fitness. The first time, the fitness measure applied was a weighted combination of the mean absolute error (MAE, weighted at 0.6) and the mean square error (MSE, weighted at 0.4). For the second training run the Structure Similarity (SSIM) index was used (Wang & Bovik 2002). The SSIM is a measure developed for use in estimating the effect of subjective viewing of structural integrity. In all other respects the training run details were identical.

	Figure 4: Output after 35 iterations using the MAE/MSE hybrid quality measure.
Open with DEXTER

The training runs were initially set for 35 iterations. Figure 4 shows the improvement, after 35 iterations, made on two examples of the training data by a filter developed using the MAE/MSE hybrid quality measure. Clearly the filter has improved these images, but there is still significant noise in one of the images. The results from the filter produced using the MAE/MSE hybrid as the error measure, with 500 iterations, is shown in Fig. 5. The final error measure at the end of this training run was 0.996. (N.B. the error measure varies from 0 to 1, with 0 being as different as possible, and 1 being a pefect match.)

	Figure 5: As in Fig. 4 but using 500 interations.
Open with DEXTER

For comparison, the results of the filter created by using the SSIM as the fitness measure is shown in Fig. 6. The final error was 0.991 and was as achieved after 500 iterations. (It is not possible to directly compare the values of the error measure obtained by different techniques.)

	Figure 6: The trained filter obtained after 500 interations using the Structure Similarity index as a quality measure.
Open with DEXTER

The filter windows determined to creates these images are shown in Figs. 7 and 8. As can be seen these two filters are very different in their operation.

	Figure 7: The hard centre and soft boundary produced by the MAE/MSE training sequence. The rank for the filter was 17 and the operation sequence was a single erosion.
Open with DEXTER

	Figure 8: The hard centre and soft boundary produced by the SSIM training model. The rank was determined to be 4, and the operation sequence was 2 erosions in sequence.
Open with DEXTER

   
5 Application to real images

The near-optimal filters produced from the two training runs were then applied to the real noisy solar images. The results using the filter trained in the first run of the application are shown in Fig. 9. The subjective improvement is very obvious. Particularly the lower-left image, which had been almost completely obscured by noise, has had a significant amount of the disturbance removed enabling solar structures to be seen. The results of applying the filter created in the second training run are shown in Fig. 10. In comparison, the second run appears to have created a better filter and has removed almost all the speckle. However there are areas in the image (the black specks) which are clearly artefacts which have been introduced, possibly due to damaged CCD pixels. Nonetheless the image is otherwise of a very high subjective quality.

	Figure 9: Two images from the LASCO C2 telescope, and their processed counterparts. This shows the result of applying the MAE/MSE combination quality measure. Some "speckle'' remains in the processed images.
Open with DEXTER

	Figure 10: The same two images as in Fig. 9 and their processed counterparts, processed using the SSIM quality measure. The "speckle'' that was present in Fig. 9 has disappeared, but some artefacts (black specks) have been introduced.
Open with DEXTER

Without clean versions of the image, with identical underlying data, measures such as the MAE, MSE and SSIM cannot be used to assess the objective improvement.

For comparison, in Figs. 11 and 12 we show the results obtained using some standard processing and cosmic-ray removal software currently available in the IDL and SolarSoft libraries. Of these, the frame-comparison method is the most successful although it is unable to significantly improve the very noisy image (Fig. 12). It is readily seen that the soft morphological processing gives qualitatively superior results to any of the other available methods.

Figure 11: The image used in the upper panel of Fig. 9 and its counterparts processed using standard techniques currently available to the community. From left to right along the top row are (1) raw image; (2) using the IDL morphological closing routine morph_close with a $[3\times 3]$ "cross of St. George'' kernel; (3) using SolarSoft zunspike.pro with the default for LASCO images. From left to right along the bottom row are (4) using median-filtering; (5) using the SolarSoft despiking algorithm nospike.pro by B. DePontieu with the same kernel as in (2), and; (6) using M. Aschwanden's SolarSoft routine trace_unspike_time.pro routine, employing interframe comparisons.

Open with DEXTER

Figure 12: The image used in the lower panel of Fig. 9 and its counterparts processed using standard techniques currently available to the community. From left to right along the top row are (1) raw image; (2) using the IDL morphological closing routine morph_close with a $[3\times 3]$ "cross of St. George'' kernel; (3) using SolarSoft zunspike.pro with the default for LASCO images. From left to right along the bottom row are (4) using median-filtering; (5) using the SolarSoft despiking algorithm nospike.pro by B. DePontieu with the same kernel as in (2), and; (6) using M. Aschwanden's SolarSoft routine trace_unspike_time.pro routine, employing interframe comparisons.

Open with DEXTER

  
6 Discussion and conclusions

This paper has presented an overview on how in practice to design and apply a soft morphological filter for the removal of cosmic ray "noise'' from astronomical images. The design and use of a training algorithm is key to the creation of a successful non-linear filter, but the time investment in this initial process is repaid by significantly improved images. In the case of the LASCO images used as examples in this paper, the improvement during epochs of heavy particle bombardment is particularly important, as it permits the morphological changes in the corona to be examined even during the earliest phases of the event that accelerated them.

Soft morphological filters deal very well with both saturated and non-saturated noise spikes, as shape as well as intensity information is incorporated. In contrast to linear operators the filters are based on a rank ordering operation with a single value being selected to form the filter output. Therefore outliers such as noise spikes are forced to the extremes of the rank ordering and have almost no chance of being selected as the output value. For linear operators, the contribution from large noise spikes would always have some weighted effect on the output. However, for these non-linear operators even large spikes are completely filtered out and do not influence the output at all.

We have tested the process on LASCO coronagraph data, and found it to function well, better than any of the existing methods we are aware of. However we note that in the case of spike removal from LASCO data, the morphological properties of the features of interest (generally smooth, linear, faint and and large-scale) and the noise (sharp, bright and small-scale) are rather different, and do not provide the most challenging test of the method. Therefore in future work we will be investigating the design and application of soft morphological filters to solar image data with smaller intrinsic scales, such as TRACE or forthcoming Solar-B data. For many applications in solar physics, the recovery of previously-obscured structures is the most important scientific output, but in other applications the absolute values of the counts per pixel are of prime importance. Therefore, further work is also necessary to quantify the photometric accuracy of the images resulting from application of the soft morphological filter.

Acknowledgements

We would like to thank our referee, whose insightful and useful comments have improved this paper. L.F. acknowledges helpful discussion with Josef Khan. The SOHO/LASCO data used here are produced by a consortium of the Naval Research Laboratory (USA), Max-Planck-Institut für Aeronomie (Germany), Laboratoire d'Astronomie (France), and the University of Birmingham (UK). We are grateful to the SOHO/LASCO team for making their data publicly available. SOHO is a project of international cooperation between ESA and NASA.

References

• Anderson, T. W. 1993, An Introduction to Multivariate Statistical Analysis, 3rd edition (New York: Wiley)
• Angulo, J., & Serra, J. 2003, in ICIP, 2003, 125 (In the text)
• Brueckner, G. E., Howard, R. A., Koomen, M. J., et al. 1995, Sol. Phys., 162, 357 [NASA ADS] [CrossRef] (In the text)
• Delaboudinière, J.-P., Artzner, G. E., Brunaud, J., et al. 1995, Sol. Phys., 162, 291 [NASA ADS] [CrossRef] (In the text)
• Domingo, V., Fleck, B., & Poland, A. I. 1995, Sol. Phys., 162, 1 [NASA ADS] [CrossRef] (In the text)
• Farrage, C. L., & Pimblett, K. A. 2005, Proc. Astron. Soc. Austr., 22, 249 [NASA ADS]
• Hamid, M. S., Marshall, S., & Harvey, N. 2003, IEEE Trans. Circuits and Systems for Video Technology, 13, 406 (In the text)
• Handy, B. N., Acton, L. W., Kankelborg, C. C., et al. 1999, Sol. Phys., 187, 229 [NASA ADS] [CrossRef] (In the text)
• Harvey, N., & Marshall, S. 1996, in Mathematical Morphology and its Applications to Image and Signal Processing, ed. P. Maragos, R. W. Schafer, & M. A. Butt (Kluwer), 186
• Harvey, N., & Marshall, S. 2000, in Mathematical Morphology and its Applications to Image Analysis, ISMM 2000, ed. J. Goutsias, L. Vincent, & D. S. Bloomberg (Kluwer) 129 (In the text)
• Heijmans, H. J. 1994, Morphological Operators (New York: Academic Press) (In the text)
• Holland, J. H. 1995, Adaptation in Natural and Artificial Systems (MIT Press), 89 (In the text)
• Kahler, S. W. 1992, ARA&A, 30, 113 [NASA ADS] (In the text)
• Koskinen, L., & Astola, J. 1994, J. Electron. Imag., 3, 60 [NASA ADS] (In the text)
• Kraft, P., Harvey, N., & Marshall, S. 1997, J. Electron. Imag., 6, 504 [NASA ADS] (In the text)
• Marshall, S., Harvey, N., & Greenhalgh, D. 2000, in X International Signal Processing Conference (EUSIPCO 2000), 392 (In the text)
• Matheron, G. 1975, Random Sets and Integral Geometry (New York: Wiley) (In the text)
• Ogawara, Y., Takano, T., Kato, T., et al. 1991, Sol. Phys., 136, 10 (In the text)
• Shaw, R. A., & Horne, K. 1992, Astron. Soc. Pacific, 25, 311 [NASA ADS] (In the text)
• Serra, J. 1982, Image Analysis and Mathematical Morphology (London: Academic Press) (In the text)
• Serra, J. 1988, Image Analysis and Mathematical Morphology, 2 (London: Academic Press)
• Wang, Z., & Bovik, A. C. 2002, IEEE Sign. Proc. Lett., 9, 81 (In the text)