An astronomical imaging stripe noise correction method and system based on low-frequency directivity features

Abstract

The application discloses a kind of astronomical imaging stripe noise correction method and system based on low-frequency directionality feature, first, realize column level brightness equalization through robust statistics, preliminarily eliminate global vertical stripe;Again, through the frequency domain analysis of sliding window, extract low-frequency directionality feature and generate heat map, obtain fixed trapezoidal mask after threshold segmentation and morphological processing;Finally, under the constraint of mask, through star point suppression, frequency domain filtering generates stripe template, and the image is compensated and corrected.This method combines spatial equalization and frequency domain template estimation, and accurately constrains the effective correction area through trapezoidal mask, which not only retains the true star point signal, but also efficiently suppresses residual stripe noise, significantly improves the signal-to-noise ratio of astronomical image and the reliability of subsequent celestial body detection.

Description

Technical Field

[0001] This invention relates to the field of astronomical survey data processing and image quality enhancement technology, specifically to an astronomical imaging stripe noise correction method and system based on low-frequency directional characteristics. It is applicable to the automated preprocessing of multi-channel stitched CCD imaging data, and is especially applicable to scientific images (SCI) acquired by sky survey cameras. Background Technology

[0002] Astronomical sky surveys continuously acquire large-scale sky image data through wide-field imaging systems, providing a crucial data foundation for celestial object detection, photometry and positioning, transient source discovery, and subsequent physical parameter inversion. In the data processing chain of survey cameras, scientific images (SCI) typically refer to the raw observation images (containing celestial signals and background noise) obtained directly after exposure. Their quality directly impacts the performance limits of key tasks such as source detection, background estimation, image overlay, and differencing. Therefore, image preprocessing and noise correction for SCI are crucial prerequisites for ensuring the scientific output of sky survey data.

[0003] As detector size and readout rates increase, sky survey cameras are increasingly adopting multi-channel stitched CCD architectures to achieve larger fields of view and higher sampling efficiency. These systems typically consist of multiple CCD devices stitched together to form a large array, employing multi-channel parallel readout to meet high-throughput observation requirements. In engineering practice, it is difficult to achieve perfect consistency in the gain, bias, and noise characteristics of different readout channels. Furthermore, factors such as electronic crosstalk, clock interference, power supply ripple, and temperature drift can affect the SCI (Search Channel Interference), leading to structured noise with a clear directionality. One typical manifestation of this is stripe noise extending along the row / column direction.

[0004] Stripe noise is characterized by its "clear directionality, prominent low-frequency components, and non-stationary spatial distribution." On the one hand, stripes can manifest as overall brightness shifts between columns or slowly undulating superimposed local periodic structures; on the other hand, their intensity and shape may vary with position, becoming more pronounced in areas close to boundaries or in certain channel regions. This type of structured noise significantly interferes with background statistical characteristics, causing an increase in background RMS and enhanced noise correlation, which in turn leads to problems such as unstable threshold detection and photometric zero-point shift, especially in scenes with weak signals, complex background undulations, or dense star points.

[0005] In existing survey data processing workflows, conventional calibration and correction (such as baseline offset correction) are typically performed first to reduce the impact of readout baseline drift. However, directional fringe structures may still remain in actual survey SCI. Common engineering processing strategies for fringe noise can be broadly categorized into two types: one is column-level correction based on statistical alignment, which estimates the background representative value of each column (or narrow column band) and aligns it to a unified reference to reduce overall inter-column offset; the other is fringe processing based on frequency domain analysis, which locates the directional components corresponding to the fringes in the spectrum and performs filtering or template estimation, then returns to the image domain for image compensation and correction. While these methods can improve image quality in specific scenarios, they still have limitations under the complex conditions of survey SCI. For example, column-level statistics are easily affected by bright stars, dense star points, and abnormal bright spots, resulting in unstable estimation; uniform frequency domain filtering across the entire image is difficult to balance the positioning accuracy of local fringes with background fidelity; at the same time, the bright star point structure generates strong energy components in the frequency domain, and if there is a lack of effective regional constraints and interference suppression mechanisms, the fringe template estimation is easily contaminated, leading to under-correction, over-correction, or the introduction of new artifacts.

[0006] In summary, while existing methods can reduce the impact of stripe noise to some extent, they still have shortcomings in dealing with common scenarios in multi-channel stitched CCD survey SCI, such as local non-stationary stripes, close boundaries, background undulations, and star point high-brightness interference. Summary of the Invention

[0007] To address the common problems of directional fringe noise and low-frequency background undulations in scientific images (SCI) from sky survey cameras, this invention proposes a method and system for astronomical imaging fringe noise correction based on low-frequency directional features. By introducing key steps such as column-level robust brightness alignment, sliding window low-frequency directional feature thermal map localization, fixed trapezoidal geometric mask constraints, and frequency domain fringe template compensation correction, the method improves the stability, controllability, and batch processing efficiency of fringe correction while ensuring that the true structural information of celestial objects is not affected as much as possible. This provides reliable technical support for high-quality preprocessing and subsequent scientific measurements of large-scale sky survey SCI.

[0008] To achieve the above-mentioned technical objectives, the technical solution adopted by the present invention is as follows:

[0009] In a first aspect, the present invention discloses an astronomical imaging stripe noise correction method based on low-frequency directional characteristics, the method comprising the following steps:

[0010] S1. Acquire two-dimensional astronomical image data collected by the sky survey camera, use robust statistical estimation to estimate the representative brightness of each column and perform column-level brightness translation, perform vertical stripe brightness equalization on the two-dimensional astronomical image, and obtain the two-dimensional astronomical image after column-level brightness equalization.

[0011] S2, with a set window size and step size, slides and samples on the two-dimensional astronomical image after column-level brightness equalization. For each window, it performs mean removal, windowing, and two-dimensional Fourier transform to extract low-frequency ring energy and statistically analyzes the directional energy distribution by angle binning. It calculates low-frequency directivity coefficients and weights the results of each window according to the coverage area to generate a low-frequency directivity feature heatmap of the same size as the original image to characterize the spatial position and relative intensity of residual fringes.

[0012] S3, after determining the contrast validity of the low-frequency directional feature heatmap, perform threshold segmentation, morphological closing operation and denoising, and retain only the connected regions with the largest area and that meet the area ratio constraint as the main region mask; fit the left and right boundary lines of the main region mask line by line, and intersect the two lines with the upper and lower boundaries of the image to obtain four vertices, generating a fixed trapezoidal mask.

[0013] S4 performs star point suppression or star point masking on the original two-dimensional astronomical image, then performs a two-dimensional Fourier transform on the processed two-dimensional astronomical image and retains only the vertical frequency band near the center of the frequency domain to obtain vertical stripes / low-frequency components. The inverse transform is then used to obtain a two-dimensional stripe template in the image domain. Under the constraint of a fixed trapezoidal mask, the two-dimensional stripe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation and correction, and the image after stripe noise correction is output.

[0014] Furthermore, in step S1, the process of performing vertical stripe brightness equalization on the two-dimensional astronomical image includes the following steps:

[0015] Divide the two-dimensional astronomical image into multiple vertical bars along the column direction, and then... A fixed number of sample points are randomly selected from each vertical bar without replacement to form a sample set. The median of the sample set is first calculated as the robust central value. Then calculate the absolute median and convert it to the Gaussian equivalent standard deviation. Only those falling within the main peak area are retained. The sample within is used as the first The set of main peaks in the vertical bars, parameters Used to remove outliers, including star-bright pixels, bad pixels, or extreme noise.

[0016] The average value of the main peak set in each vertical column is taken as the main peak brightness estimate. The main peak sequence is obtained by combining the main peak brightness estimates of all columns. The global mean of the brightness of all vertical main peaks was calculated. , Indicates the first Estimated brightness values ​​of the main peaks of the vertical bars;

[0017] global mean As the alignment target for each vertical bar, the first The correction offset for each vertical bar is:

[0018] ;

[0019] in To correct the strength coefficient, Indicates perfect alignment. Partial correction is performed to reduce the risk of overcorrection; then, based on the correction offset, brightness translation correction is performed on all pixels in each vertical bar to obtain a brightness-equalized image, thus obtaining a two-dimensional astronomical image with column-level brightness equalization.

[0020] Furthermore, in step S2, the process of generating a low-frequency directional feature heatmap includes the following steps:

[0021] Let the size of the sliding window be... Step size is , and These represent the number of pixels the sliding window moves vertically and horizontally each time; a sequence of starting points at the top left corner of the window is generated using an edge-fitting overlay strategy, and for any window starting point... Take a sub-block The window is subjected to mean removal and windowing to obtain preprocessed sub-blocks. :

[0022] ;

[0023] in The mean of the window. This represents element-wise multiplication. A two-dimensional Hann window; based on preprocessed sub-blocks Calculate the frequency domain result after two-dimensional Fourier transform And the power spectrum was calculated. In the formula, Represents a two-dimensional Fourier transform operator; This represents a two-dimensional astronomical image after column-level brightness equalization.

[0024] The power spectrum is centered, and the radius is defined with the center of the spectrum as the origin. Select the low-frequency ring zone region As a statistical object, and Evenly divided into Each angle is divided into boxes, and the energy within the ring is aggregated according to the angle to obtain directional energy. The low-frequency directivity coefficient is calculated using the following formula:

[0025] ;

[0026] in To prevent the stable term from being divided by zero; The directional energy obtained by aggregation within the k-th angle sub-bin;

[0027] Initialize the cumulative graph With the hit chart For each window, the low-frequency directional coefficients are set to zero. Backfill to all pixel locations within the window's coverage area, so that After all windows have been processed, normalization is performed on each pixel to obtain a low-frequency directional feature heatmap H:

[0028] ;

[0029] In the formula, Indicates pixel position The normalized value of the low-frequency directional characteristic heatmap at that location. Indicates pixel position Above, the low-frequency directivity coefficients of all sliding windows that cover it. The sum, Indicates the location of the overlaid pixel. The total number of sliding windows in the area.

[0030] Furthermore, in step S2, multiple window patches are batch-assembled into tensors and FFT and directional energy aggregation are performed in parallel on the GPU.

[0031] Step S3 further includes:

[0032] Calculate the low-frequency directional characteristic heatmap dynamic range with standard deviation :

[0033] ;

[0034] In the formula, and These represent heatmaps showing low-frequency directional characteristics. The minimum and maximum values ​​of all pixels in the array. This represents the function for calculating standard deviation.

[0035] Heatmap for determining low-frequency directional characteristics Whether the dynamic range and standard deviation are less than their respective preset thresholds, when or hour, The dynamic range threshold and standard deviation threshold are used respectively. If the heatmap lacks effective contrast, an empty mask is directly returned. Otherwise, if the low-frequency directional feature heatmap meets the contrast requirement, it is considered to be smoothed to suppress local noise fluctuations. Then, threshold segmentation is used to verify the candidate regions, resulting in candidate binary maps. :

[0036] ;

[0037] In the formula, Indicates the threshold; Indicates pixel position Binary candidate region value at;

[0038] For candidate binary graphs After performing morphological closing and hole-filling operations, small connected components are removed, connected component labeling is performed, and only the connected component with the largest area is retained as the main region mask. Calculate the main region mask area ratio :

[0039] ;

[0040] when or Return to empty mask at time, and These represent the height and width of a two-dimensional astronomical image, respectively. and These represent the minimum and maximum thresholds for the area proportion of the main region, respectively. Indicates the main region mask The number of pixels contained; when returning an empty mask, directly output the two-dimensional astronomical image after column-level brightness equalization, and end the process.

[0041] Furthermore, the threshold Use Otsu adaptive thresholding or percentile thresholding to obtain the threshold.

[0042] Furthermore, in step S3, the process of generating the fixed trapezoidal mask includes the following steps:

[0043] masking of the main region Perform a fixed trapezoidal fit. Specifically, search the set of pixel column indices where the mask is true for each row, and take the leftmost and rightmost indices as the left and right boundary observation points. and , Represents the row coordinates of the image. and They represent the first In-line main region mask The leftmost and rightmost column coordinates of the true pixels are given; a linear fit is then performed on the left and right boundary observation points respectively.

[0044] ;

[0045] In the formula, and These represent the slope and intercept of the straight line fitted from the observation points on the left boundary, respectively. and These are the slope and intercept of the straight line fitted from the observation points on the right boundary, respectively.

[0046] Then, force the ordinate of the trapezoid's vertices to be set to the upper and lower boundaries of the image. and This gives us the x-coordinates of the four vertices.

[0047] ;

[0048] And cut it to The legal range, in the formula, , , and Let x and y represent the x-coordinates of the top-left, bottom-left, top-right, and bottom-right vertices of the trapezoid, respectively. and These represent the height and width of a two-dimensional astronomical image, respectively. When the x-coordinate of a vertex is 0 or... Less than hour, As the adsorption ratio parameter, this vertex is attached to 0 or... Finally, a trapezoidal polygon is constructed using the four vertices, and a fixed trapezoidal mask is generated. .

[0049] Furthermore, in step S4, the process of obtaining the two-dimensional stripe template in the image domain includes the following steps:

[0050] For a two-dimensional astronomical image after column-level brightness equalization, the set of pixels in a certain direction Computational Robustness Center

[0051] ;

[0052] and MAD equivalent standard deviation

[0053] ;

[0054] In the formula, This represents the function for calculating the median.

[0055] When the pixel value of any pixel When this happens, the pixel is recorded as a star point and written into the mask, while the pixel value is replaced with the median or neighborhood statistics in the corresponding direction. The threshold coefficient for star point discrimination is represented; after iteration and morphological dilation, the mask region is filled with the neighborhood median or edge statistics to obtain a spatially continuous star point suppression background map. ;

[0056] Suppressing background image with star points The global mean is removed to suppress the DC component, and after multiplying by a two-dimensional Hann window, a two-dimensional Fourier transform is calculated and centered to obtain the centered frequency domain image. :

[0057] ;

[0058] In the formula, Represents the two-dimensional Fourier transform operator. This represents the Fourier transform centering operator;

[0059] Let the frequency domain center column be , Represents the width of a two-dimensional astronomical image, preserving the column range. ,in The half-width parameter is used to determine the frequency range covered by the template and to construct the frequency domain mask. Make

[0060] ;

[0061] Frequency domain preservation is then performed to obtain the frequency domain image after frequency domain masking. :

[0062] ;

[0063] By performing inverse centering and inverse Fourier transform on it, a two-dimensional stripe template in the image domain can be obtained:

[0064] ;

[0065] in The operator for taking the real part is... This represents the two-dimensional inverse Fourier transform operator. This represents element-wise multiplication; These represent the row and column coordinates of the frequency domain image, respectively. Indicates the frequency domain mask in coordinates Pixel value at that location This represents the inverse centering operator of the Fourier transform.

[0066] Furthermore, in step S4, under the constraint of a fixed trapezoidal mask, the two-dimensional stripe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation and correction, and the image after stripe noise correction is output:

[0067] ;

[0068] in The template compensation strength coefficient is used to adjust the correction amplitude. To fix the trapezoidal mask, Represents a two-dimensional stripe template in the image domain. This represents a two-dimensional astronomical image after column-level brightness equalization.

[0069] Secondly, the present invention discloses an astronomical imaging stripe noise correction system based on low-frequency directional characteristics, the system comprising:

[0070] The vertical stripe brightness equalization module is used to perform vertical stripe brightness equalization on the acquired two-dimensional astronomical image data by using robust statistical estimation of the representative brightness of each column and performing column-level brightness shifting. This results in a column-level brightness equalized two-dimensional astronomical image.

[0071] The low-frequency directional feature heatmap calculation module is used to slide sample on the two-dimensional astronomical image after column-level brightness equalization with a set window size and step size. For each window, it performs mean removal, windowing and two-dimensional Fourier transform to extract low-frequency ring energy and statistically analyze the directional energy distribution by angle binning. It calculates the low-frequency directional coefficient and weights the results of each window according to the coverage area to generate a low-frequency directional feature heatmap of the same size as the original image to characterize the spatial position and relative intensity of residual fringes.

[0072] The trapezoidal fitting module is used to determine the contrast validity of the low-frequency directional feature heatmap and then perform threshold segmentation, morphological closing operation and denoising. Only the connected regions with the largest area and that meet the area ratio constraint are retained as the main region mask. The left and right boundary lines are fitted to the main region mask line by line, and the intersection of the two lines with the upper and lower boundaries of the image is used to obtain four vertices, generating a fixed trapezoidal mask.

[0073] The stripe template estimation and compensation correction module is used to perform star point suppression or star point masking on the original two-dimensional astronomical image, and then perform a two-dimensional Fourier transform on the processed two-dimensional astronomical image to retain only the vertical frequency band near the center of the frequency domain to obtain the vertical stripes / low-frequency components. The inverse transform is used to obtain the two-dimensional stripe template in the image domain. Under the constraint of a fixed trapezoidal mask, the two-dimensional stripe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation correction, and the stripe noise corrected image is output.

[0074] Compared with the prior art, the beneficial effects of the present invention are as follows:

[0075] The astronomical imaging fringe noise correction method and system based on low-frequency directional features of this invention adopts a processing flow of "column-level brightness alignment - low-frequency directional positioning - mask constraint - fringe template estimation - template compensation correction". It can effectively reduce the interference of the star point high brightness structure on the fringe template estimation without changing the original star point signal organization form, and improve the stability and controllability of the template estimation. At the same time, by leveraging the spatial positioning capability of the low-frequency directional feature heatmap and the geometric consistency constraint of trapezoidal fitting, it enhances the robustness to scenarios with close boundaries, scale changes and background undulations, and can improve processing efficiency under batch Fourier calculation (supporting GPU acceleration). Thus, it provides a reliable technical means for high-quality preprocessing and subsequent scientific measurement of large-scale sky survey images. Attached Figure Description

[0076] Figure 1 This is a schematic diagram of the overall process of the astronomical imaging stripe noise correction method based on low-frequency directional characteristics of the present invention.

[0077] Figure 2 This is a schematic diagram comparing the brightness of the vertical stripes before and after equalization; where (a) corresponds to the original HDU1 image, (b) corresponds to the corrected HDU1 image, (c) corresponds to the original HDU2 image, and (d) corresponds to the corrected HDU2 image.

[0078] Figure 3 This is a schematic diagram showing the change of low-frequency directional features with the position of the sliding window, used to illustrate the change curve of low-frequency fringe intensity with different sliding window positions; where (a) corresponds to image example one and (b) corresponds to image example two.

[0079] Figure 4 A schematic diagram of main region extraction and fixed trapezoidal mask fitting; where (a) is the thermal... Figure 1 (b) represents heat. Figure 1 The mask, (c) is thermal Figure 2 (d) represents heat. Figure 2 The mask;

[0080] Figure 5 This diagram illustrates the results of stripe noise template extraction and correction, used to show the comparison between the generation of the stripe template and the effect of compensation and correction on the original image. (a) corresponds to the original image, (b) corresponds to the star mask, (c) corresponds to the image after 2D Fourier transform, (d) corresponds to the image after retaining the low-frequency signal, (e) corresponds to the noise template extracted by inverse Fourier transform, and (f) corresponds to the image after subtracting the noise template from the original image. Detailed Implementation

[0081] The embodiments of the present invention will be described in further detail below with reference to the accompanying drawings.

[0082] This invention discloses a method for correcting astronomical imaging stripe noise based on low-frequency directional characteristics, the method comprising the following steps:

[0083] S1. Acquire two-dimensional astronomical image data collected by the sky survey camera, use robust statistical estimation to estimate the representative brightness of each column and perform column-level brightness translation, perform vertical stripe brightness equalization on the two-dimensional astronomical image, and obtain the two-dimensional astronomical image after column-level brightness equalization.

[0084] S2, with a set window size and step size, slides and samples on the two-dimensional astronomical image after column-level brightness equalization. For each window, it performs mean removal, windowing, and two-dimensional Fourier transform to extract low-frequency ring energy and statistically analyzes the directional energy distribution by angle binning. It calculates low-frequency directivity coefficients and weights the results of each window according to the coverage area to generate a low-frequency directivity feature heatmap of the same size as the original image to characterize the spatial position and relative intensity of residual fringes.

[0085] S3, after determining the contrast validity of the low-frequency directional feature heatmap, perform threshold segmentation, morphological closing operation and denoising, and retain only the connected regions with the largest area and that meet the area ratio constraint as the main region mask; fit the left and right boundary lines of the main region mask line by line, and intersect the two lines with the upper and lower boundaries of the image to obtain four vertices, generating a fixed trapezoidal mask.

[0086] S4 performs star point suppression or star point masking on the original two-dimensional astronomical image, then performs a two-dimensional Fourier transform on the processed two-dimensional astronomical image and retains only the vertical frequency band near the center of the frequency domain to obtain vertical stripes / low-frequency components. The inverse transform is then used to obtain a two-dimensional stripe template in the image domain. Under the constraint of a fixed trapezoidal mask, the two-dimensional stripe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation and correction, and the image after stripe noise correction is output.

[0087] This invention uses a two-dimensional image array as input and employs a collaborative process of "column-level robust brightness alignment - low-frequency directional feature heatmap localization - trapezoidal geometric constraints - frequency domain vertical template estimation - template compensation correction" to achieve stable correction of directional stripes and low-frequency structures.

[0088] In the vertical stripe brightness equalization stage, this invention corrects the brightness shift between columns through column-level robust main peak estimation and a global alignment strategy. The core idea is to divide the image into columns, using sampling and robust statistics to suppress the influence of outliers such as star highlights and bad pixels on the estimation, obtaining the representative background main peak brightness for each column. Then, based on the global reference brightness, the column-level correction offset is calculated, and brightness translation correction is performed on the entire column. This weakens the baseline component of the vertical stripes without excessively damaging the true structural information, providing a more stable input for subsequent heatmap localization and frequency domain template estimation.

[0089] Specifically, let the input SCI image be a two-dimensional matrix. ,in For image height, The image width is given by pixel coordinates. To ensure the accuracy of numerical calculations, the original image is first converted to a floating-point type. Then, vertical stripe brightness equalization is performed to correct the overall brightness shift between columns. Specifically, the image is divided into columns (or narrow stripes), and the first... The set of column pixels is denoted as To reduce the computational burden of traversing the entire column and avoid a small number of abnormal bright spots dominating the statistics, the first... The column randomly selects a fixed number of sample points from the pixel set without replacement to form the sample set. ,in The parameter `samples_per_strip` specifies the robust center value. The median is first calculated for the sample set.

[0090] ;

[0091] Calculate the median absolute deviation (MAD) and convert it to the Gaussian equivalent standard deviation.

[0092] ;

[0093] in, This represents the robust center estimate of the background peak. This represents a dispersion estimate that is insensitive to outliers. This represents the grayscale value of the j-th sample pixel randomly selected from the i-th column. Then, only pixels falling within the main peak region are retained. The samples within are used as the main peak set

[0094] ;

[0095] in The parameter `gauss_k` is used to remove outliers such as bright pixels, bad pixels, or extreme noise. The average of the set of main peaks is used as the brightness estimate for that column's main peak.

[0096] ;

[0097] The system will be deactivated if the number of valid samples in the main peak set is insufficient. To ensure the stability of the estimation, the above sampling and main peak estimation process is repeated for all columns to obtain the main peak sequence. And calculate the global reference brightness (global main peak mean).

[0098] ;

[0099] Using this as the alignment target for each column, calculate the first... Column correction offset

[0100] ;

[0101] in To correct the strength coefficient, Indicates perfect alignment. This indicates partial correction to reduce the risk of overcorrection. Subsequently, brightness shift correction is performed on all pixels in the column to obtain a brightness-balanced image.

[0102] ;

[0103] In the formula, This represents the pixel value in the y-th row and i-th column of the image after brightness equalization. This represents the pixel value in the y-th row and i-th column of the floating-point image.

[0104] To prevent calibration from causing brightness to exceed limits or resulting in oversaturation, you can... The brightness range of the input image is uniformly cropped back to ensure that the output brightness range is consistent with the input. Figure 2 This diagram illustrates the comparison of vertical fringe brightness before and after column-level equalization. HDU1 and HDU2 correspond to image data of the observed target at two different bands, with different filters, or different exposures. In the low-frequency directional feature heatmap stage, this invention spatially locates residual directional fringes that may still exist after column-level alignment. The core idea is to perform a sliding window scan of the image, conduct frequency domain directional analysis within each window to obtain a window-level measure reflecting the strength of low-frequency directionality, and generate a heatmap of the same size as the original image using a fusion strategy of "window coverage backfilling—hit counting—pixel-by-pixel averaging." By configuring parameters such as window size, step size, and angle binning, a balance can be achieved between spatial positioning accuracy and result smoothness and stability, thus enabling the heatmap to highlight significant fringe areas while avoiding local misjudgments caused by random noise. This invention achieves automatic fringe region localization through low-frequency directional feature heatmaps, avoiding undercorrection or overcorrection caused by uniform strong processing of the entire image.

[0105] After completing column-level brightness equalization, further steps are taken to locate the spatial distribution of residual stripes or directional low-frequency background undulations. A sliding window frequency domain analysis is performed to generate a low-frequency directional characteristic heatmap. Let the sliding window size be... Step size is A "edge-covering" strategy is used to generate the starting sequence at the top left corner of the window, ensuring that the right and bottom boundary regions are scanned and covered, thus avoiding unanalyzed edge bands. For any window starting point... Take a sub-block

[0106] ;

[0107] To suppress the DC component and reduce spectral leakage, mean removal and windowing are performed on the window to obtain...

[0108] ;

[0109] in The mean of the window. This represents element-wise multiplication. A two-dimensional Hann window is used (obtained by the outer product of two one-dimensional Hann windows). Then, the two-dimensional Fourier transform and power spectrum are calculated.

[0110] ;

[0111] The power spectrum is centered (fftshift), and the radius is defined with the center of the spectrum as the origin. Select the low-frequency ring zone region As a statistical object, and Evenly divided into Each angle is divided into boxes, and the energy within the ring is aggregated according to the angle to obtain directional energy. Therefore, the low-frequency anisotropy coefficient (which can also be expressed as the low-frequency directivity coefficient) is defined.

[0112] ;

[0113] in To prevent the stable term from being divided by zero, A larger value indicates that low-frequency energy is more concentrated in a certain direction, usually corresponding to more obvious stripes or directional textures. To obtain a heatmap of the same size as the original image, this embodiment employs a fusion mechanism of "window overlay backfill—hit count—pixel-by-pixel averaging": initializing the cumulative image... With the hit chart For a zero matrix, obtain a scalar for each window. Then it is backfilled to all pixel locations within the window's coverage area, so that

[0114] ;

[0115] After all windows have been processed, normalization is performed on each pixel to obtain the heatmap.

[0116] ;

[0117] The more windows overlap (the smaller the step size), the better. The larger the step size, the smoother and more stable the heatmap, but the computational cost increases; a larger step size emphasizes local responses, resulting in a more pronounced blocky scale in the heatmap. To improve window-level computation efficiency, multiple window patches can be batch-assembled into tensors and FFT and directional energy aggregation can be performed in parallel on the GPU, thereby significantly improving the heatmap generation speed. Figure 3 This is a schematic diagram showing the change of low-frequency directional features with the position of the sliding window, used to illustrate the change curve of low-frequency fringe intensity with different sliding window positions; where (a) corresponds to image example one and (b) corresponds to image example two.

[0118] In the main region extraction and fixed trapezoidal fitting stages, this invention transforms the main abnormal regions in the heatmap into geometrically consistent and reusable effective processing region constraints. The core idea is to first determine the contrast validity of the heatmap to avoid unreliable region fitting in low-contrast scenes; if valid, candidate regions are obtained through threshold segmentation and morphological processing, retaining only the largest connected component as the main region; then, the left and right boundary points of the main region are extracted row by row and linearly fitted, ultimately constructing a fixed trapezoidal mask. By introducing rules such as vertex edge snapping, the geometric stability when the main region is close to the image boundary can be enhanced, avoiding unnecessary mask drift near the boundary, thereby improving the consistency and controllability of subsequent template estimation and correction. This invention provides geometric consistency constraints through a fixed trapezoidal mask, improving robustness in scenes with close boundaries, scale changes, and background undulations.

[0119] This step is mainly based on heat maps. The main anomaly region is automatically extracted and a fixed trapezoidal mask is fitted to impose geometric constraints on subsequent stripe template estimation. Specifically, to avoid threshold segmentation failure due to insufficient overall contrast of the heatmap, its dynamic range and standard deviation are first calculated.

[0120] ;

[0121] when or ( If the threshold is not met, the heatmap is considered to lack effective contrast, and an empty mask is directly returned to avoid unreliable region fitting. When the contrast requirement is met, the heatmap can be smoothed to suppress local noise fluctuations, and then threshold segmentation is used to verify candidate regions; the threshold can be determined using Otsu adaptive thresholding or percentile thresholding; for example, the threshold T is calculated using the following formula:

[0122] ;

[0123] in For percentile parameters, () represents the percentile statistical function, which subsequently yields the binary candidate region.

[0124] ;

[0125] To enhance regional connectivity and remove fragment noise, candidate binary maps are processed. After performing morphological closing and hole-filling operations, small connected components are removed, connected component labeling is performed, and only the connected component with the largest area is retained as the main region mask. To avoid the fitting becoming meaningless due to the region being too small or almost covering the entire map, a constraint is imposed on the area ratio of the main region.

[0126] ;

[0127] when or Return to the empty mask.

[0128] After obtaining the main region mask Subsequently, to construct a geometrically consistent and reusable constraint region, this embodiment performs a fixed trapezoidal fitting on the main region. Specifically, for each row... Search the set of pixel indices where the mask is true, and take the leftmost and rightmost indices as the left and right boundary observation points. and When there are enough effective rows, perform a linear fit on the left and right boundary points respectively.

[0129] ;

[0130] Then, the ordinates of the trapezoidal vertices are forcibly set to the upper and lower boundaries of the image. and This gives us the x-coordinates of the four vertices.

[0131] ;

[0132] And cut it to The legal range. To improve the geometric stability when the main region is close to the image boundary, a rule of "snagging to the left and right boundaries" in the x-direction can be introduced: when the x-coordinate distance of a vertex is 0 or... Less than ( When the adsorption ratio parameter is used, the vertex is forced to be attached to 0 or... This avoids unnecessary drift during fitting near the boundaries. Finally, a trapezoidal polygon is constructed with four vertices, and a fixed trapezoidal mask is generated. This serves as an effective region constraint for subsequent frequency domain template estimation. Figure 4 A schematic diagram of main region extraction and fixed trapezoidal mask fitting; where (a) is the thermal... Figure 1 (b) represents heat. Figure 1 The mask, (c) is thermal Figure 2 (d) represents heat. Figure 2 The mask.

[0133] In the frequency domain fringe template estimation and compensation correction stage, this invention further corrects residual directional fringes and low-frequency structures under mask constraints. The core idea is: to reduce the interference of strong energy generated by bright star structures in the frequency domain on fringe template estimation, firstly, star suppression results or star mask processing results are constructed through row / column robust discrimination; then, a two-dimensional Fourier transform is performed on the image to construct a vertical frequency band preservation region near the center of the frequency domain to extract fringe components, and an inverse transform is performed to obtain a two-dimensional fringe template in the image domain; finally, this two-dimensional fringe template is directly used to compensate and correct the original SCI image within the trapezoidal mask constraint range, thereby weakening directional fringes and low-frequency structures while preserving the true structural information of celestial bodies. This stage can control the trade-off between correction intensity and background fidelity through parameters such as bandwidth, mean removal, and windowing to adapt to changes in fringe intensity and morphology under different observation conditions. This invention, by estimating the frequency domain vertical fringe template under mask constraints and directly compensating and correcting the original star-containing SCI image, effectively reduces the interference of bright star structures on template estimation and improves the controllability and repeatability of the correction results.

[0134] To avoid strong energy interference from bright structures such as star points in the frequency domain, which could affect the stable estimation of the fringe template, a star point suppression background map is first constructed only along the "template estimation path". This image is used as a template for subsequent FFT estimation, but is not a direct replacement image for the final output. Star suppression can be implemented using row / column robustness. -Clipping method: For a set of pixels in a certain direction Computational Robustness Center

[0135] ;

[0136] and MAD equivalent standard deviation

[0137] ;

[0138] and with threshold Identify star points / abnormally bright pixels when pixel value The star points are recorded as star points and written into the mask. Simultaneously, the pixel value is replaced in the working image with the median or neighborhood statistics of the corresponding direction. This process can be iterated multiple times to enhance robustness to strong star points and complex backgrounds. Morphological dilation of the mask can be applied to cover the outer edges of the star points and bright halo areas. Subsequently, the mask area is filled with background using the mean of the inline boundary median or local median filtering, thus obtaining a spatially continuous star-point suppressed background image. It is important to emphasize that star suppression / filling is only used to improve the robustness of frequency domain template estimation. The final correction will be applied to the original star-containing image, thereby avoiding the introduction of additional errors by "back-injecting" or reconstructing the real star signal.

[0139] Subsequently The vertical fringe template is extracted in the frequency domain. Optionally, the global mean is removed first to suppress the DC component, and a two-dimensional Hann window can be applied to reduce spectral leakage caused by boundary discontinuities. Then, the two-dimensional Fourier transform is calculated and centered.

[0140] ;

[0141] To extract the vertical stripe / low-frequency structure, a vertical retention band is constructed near the center of the frequency domain: let the center column of the frequency domain be... Preserve column range ,in This is the strip-half-width parameter. Construct the frequency domain mask. Make

[0142] ;

[0143] And perform frequency domain preservation to obtain

[0144] ;

[0145] Inverse centering and inverse Fourier transform are performed on it to obtain a two-dimensional stripe template in the image domain.

[0146] ;

[0147] in This represents the real part operator. (This template...) It mainly includes the low-frequency / near-DC directional structure corresponding to the central vertical band in the frequency domain, which can be used to characterize residual fringes or gently varying directional undulations in the background; parameters Determine the frequency range included in the template. The larger the template size, the more low-frequency components it contains, resulting in stronger correction, but it may also remove some of the true background gradient. The smaller the value, the more focused the fringe-related frequencies become, and the more conservative the correction.

[0148] Finally, template compensation correction is performed to obtain the output image. This embodiment uses a method of "direct compensation of the original star-containing image using a template" to complete the correction, and applies correction only to the effective region under trapezoidal mask constraints, i.e., letting...

[0149] ;

[0150] in This is the template compensation strength coefficient, used to adjust the correction amplitude. A fixed trapezoidal mask is used. This design employs star-suppressed background images to improve stability during the template estimation stage, while the final correction stage directly applies the correction to the original star-containing image. This achieves effective correction of stripe noise and directional low-frequency structures without introducing star-suppressed errors. If the heatmap contrast effectiveness determination fails or the fixed trapezoidal mask is empty, the frequency domain template correction stage can be skipped, and the column-level brightness-equalized image can be directly output. This is to avoid unnecessary processing of the image when there is no valid stripe indication. Figure 5 This diagram illustrates the results of stripe noise template extraction and correction, used to show the comparison between the generation of the stripe template and the effect of compensation and correction on the original image. (a) corresponds to the original image, (b) corresponds to the star mask, (c) corresponds to the image after 2D Fourier transform, (d) corresponds to the image after retaining the low-frequency signal, (e) corresponds to the noise template extracted by inverse Fourier transform, and (f) corresponds to the image after subtracting the noise template from the original image.

[0151] Through the above implementation methods, this invention can stably locate significant fringe regions and estimate vertical low-frequency templates under geometric constraints in large-scale sky survey SCI data processing. Then, it compensates and corrects the original star-containing image with controllable intensity, thereby reducing the impact of directional fringes and low-frequency pseudostructures on subsequent scientific measurements such as source extraction and photometry. Simultaneously, it utilizes mechanisms such as sliding window batch FFT and geometric caching to improve overall processing efficiency, meeting the requirements of engineering batch processing. Furthermore, the correction method proposed in this invention can be combined with batch computing strategies to meet the engineering requirements of efficiency and stability for automated processing of large-scale sky survey data.

[0152] The present invention also provides an astronomical imaging stripe noise correction system based on the above method. The system can be deployed on a server containing a GPU or an astronomical data processing platform to realize batch and automated stripe noise correction processing of SCI images from sky survey cameras, providing more stable input data for subsequent source extraction, photometry and positioning and scientific analysis.

[0153] Although preferred embodiments of this application have been described, those skilled in the art, upon learning the basic inventive concept, can make other changes and modifications to these embodiments. Therefore, the appended claims are intended to be interpreted as including the preferred embodiments as well as all changes and modifications falling within the scope of this application.

[0154] Obviously, those skilled in the art can make various modifications and variations to this application without departing from the spirit and scope of this application. Therefore, if such modifications and variations fall within the scope of the claims of this application and their equivalents, this application also intends to include such modifications and variations.

Claims

1. A method for correcting astronomical imaging stripe noise based on low-frequency directional characteristics, characterized in that, The method includes the following steps: S1. Acquire two-dimensional astronomical image data collected by the sky survey camera, use robust statistical estimation to estimate the representative brightness of each column and perform column-level brightness translation, perform vertical stripe brightness equalization on the two-dimensional astronomical image, and obtain the two-dimensional astronomical image after column-level brightness equalization. S2, with a set window size and step size, slides and samples on the two-dimensional astronomical image after column-level brightness equalization. For each window, it performs mean removal, windowing, and two-dimensional Fourier transform to extract low-frequency ring energy and statistically analyzes the directional energy distribution by angle binning. It calculates low-frequency directivity coefficients and weights the results of each window according to the coverage area to generate a low-frequency directivity feature heatmap of the same size as the original image to characterize the spatial position and relative intensity of residual fringes. S3, after determining the contrast validity of the low-frequency directional feature heatmap, perform threshold segmentation, morphological closing operation and denoising, and retain only the connected regions with the largest area and that meet the area ratio constraint as the main region mask; fit the left and right boundary lines of the main region mask line by line, and intersect the two lines with the upper and lower boundaries of the image to obtain four vertices, generating a fixed trapezoidal mask. S4. Perform star point suppression or star point masking on the two-dimensional astronomical image after column-level brightness equalization. Then perform a two-dimensional Fourier transform on the processed two-dimensional astronomical image and retain only the vertical frequency band near the center of the frequency domain to obtain the vertical stripes / low-frequency components. The inverse transform is used to obtain the two-dimensional stripe template in the image domain. Under the constraint of a fixed trapezoidal mask, the two-dimensional stripe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation and correction, and the image after stripe noise correction is output.

2. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 1, characterized in that, Step S1, the process of performing vertical stripe brightness equalization on the two-dimensional astronomical image, includes the following steps: Divide the two-dimensional astronomical image into multiple vertical bars along the column direction, and then... A fixed number of sample points are randomly selected from each vertical bar without replacement to form a sample set. The median of the sample set is first calculated as the robust central value. Then calculate the absolute median and convert it to the Gaussian equivalent standard deviation. Only those falling within the main peak area are retained. The sample within is used as the first The set of main peaks in the vertical bars, parameters Used to remove outliers, including star-bright pixels, bad pixels, or extreme noise. The average value of the main peak set in each vertical column is taken as the main peak brightness estimate. The main peak sequence is obtained by combining the main peak brightness estimates of all columns. The global mean of the brightness of all vertical main peaks was calculated. , Indicates the first Estimated brightness values ​​of the main peaks of the vertical bars; global mean As the alignment target for each vertical bar, the first The correction offset for each vertical bar is: ； in To correct the strength coefficient, Indicates perfect alignment. Partial correction is performed to reduce the risk of overcorrection; then, based on the correction offset, brightness translation correction is performed on all pixels in each vertical bar to obtain a brightness-equalized image, thus obtaining a two-dimensional astronomical image with column-level brightness equalization.

3. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 1, characterized in that, In step S2, the process of generating a low-frequency directional feature heatmap includes the following steps: Let the size of the sliding window be... Step size is , and These represent the number of pixels the sliding window moves vertically and horizontally each time; a sequence of starting points at the top left corner of the window is generated using an edge-fitting overlay strategy, and for any window starting point... Take a sub-block The window is subjected to mean removal and windowing to obtain preprocessed sub-blocks. : ； in The mean of the window. This represents element-wise multiplication. A two-dimensional Hann window; based on preprocessed sub-blocks Calculate the frequency domain result after two-dimensional Fourier transform And the power spectrum was calculated. In the formula, Represents a two-dimensional Fourier transform operator; This represents a two-dimensional astronomical image after column-level brightness equalization. The power spectrum is centered, and the radius is defined with the center of the spectrum as the origin. Select the low-frequency ring zone region As a statistical object, and Evenly divided into Each angle is divided into boxes, and the energy within the ring is aggregated according to the angle to obtain directional energy. The low-frequency directivity coefficient is calculated using the following formula: ； in To prevent the stable term from being divided by zero; The directional energy obtained by aggregation within the k-th angle sub-bin; Initialize the cumulative graph With the hit chart For each window, the low-frequency directional coefficients are set to zero. Backfill to all pixel locations within the window's coverage area, so that After all windows have been processed, normalization is performed on each pixel to obtain a low-frequency directional feature heatmap H: ； In the formula, Indicates pixel position The normalized value of the low-frequency directional characteristic heatmap at that location. Indicates pixel position Above, the low-frequency directivity coefficients of all sliding windows that cover it. The sum, Indicates the location of the overlaid pixel. The total number of sliding windows in the area.

4. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 3, characterized in that, In step S2, multiple window patches are batch assembled into tensors and FFT and directional energy aggregation are performed in parallel on the GPU.

5. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 1, characterized in that, Step S3 further includes: Calculate the low-frequency directional characteristic heatmap dynamic range with standard deviation : ； In the formula, and These represent heatmaps showing low-frequency directional characteristics. The minimum and maximum values ​​of all pixels in the array. This represents the function for calculating standard deviation. Heatmap for determining low-frequency directional characteristics Whether the dynamic range and standard deviation are less than their respective preset thresholds, when or hour, The dynamic range threshold and standard deviation threshold are used respectively. If the heatmap lacks effective contrast, an empty mask is directly returned. Otherwise, if the low-frequency directional feature heatmap meets the contrast requirement, it is considered to be smoothed to suppress local noise fluctuations. Then, threshold segmentation is used to verify the candidate regions, resulting in candidate binary maps. : ； In the formula, Indicates the threshold; Indicates pixel position Binary candidate region value at; For candidate binary graphs After performing morphological closing and hole-filling operations, small connected components are removed, connected component labeling is performed, and only the connected component with the largest area is retained as the main region mask. Calculate the main region mask area ratio : ； when or Return to empty mask at time, and These represent the height and width of a two-dimensional astronomical image, respectively. and These represent the minimum and maximum thresholds for the area proportion of the main region, respectively. Indicates the main region mask The number of pixels contained; when returning an empty mask, directly output the two-dimensional astronomical image after column-level brightness equalization, and end the process.

6. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 5, characterized in that, The threshold Use Otsu adaptive thresholding or percentile thresholding to obtain the threshold.

7. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 1, characterized in that, In step S3, the process of generating the fixed trapezoidal mask includes the following steps: masking of the main region Perform a fixed trapezoidal fit. Specifically, search the set of pixel column indices where the mask is true for each row, and take the leftmost and rightmost indices as the left and right boundary observation points. and , Represents the row coordinates of the image. and They represent the first In-line main region mask The leftmost and rightmost column coordinates of the true pixels are given; a linear fit is then performed on the left and right boundary observation points respectively. ； In the formula, and These represent the slope and intercept of the straight line fitted from the observation points on the left boundary, respectively. and These are the slope and intercept of the straight line fitted from the observation points on the right boundary, respectively. Then, force the ordinate of the trapezoid's vertices to be set to the upper and lower boundaries of the image. and This gives us the x-coordinates of the four vertices. ； And cut it to The legal range, where, , , and Let x and y represent the x-coordinates of the top-left, bottom-left, top-right, and bottom-right vertices of the trapezoid, respectively. and These represent the height and width of a two-dimensional astronomical image, respectively. When the x-coordinate of a vertex is 0 or... Less than hour, As the adsorption ratio parameter, this vertex is attached to 0 or... Finally, a trapezoidal polygon is constructed using the four vertices, and a fixed trapezoidal mask is generated. .

8. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 1, characterized in that, Step S4, the process of obtaining the two-dimensional stripe template in the image domain includes the following steps: For a two-dimensional astronomical image after column-level brightness equalization, the set of pixels in a certain direction Computational Robustness Center ； and MAD equivalent standard deviation ； In the formula, This represents the function for calculating the median. When the pixel value of any pixel When this happens, the pixel is recorded as a star point and written into the mask, while the pixel value is replaced with the median or neighborhood statistics in the corresponding direction. The threshold coefficient for star point discrimination is represented; after iteration and morphological dilation, the mask region is filled with the neighborhood median or edge statistics to obtain a spatially continuous star point suppression background map. ; Suppressing background image with star points The global mean is removed to suppress the DC component, and after multiplying by a two-dimensional Hann window, a two-dimensional Fourier transform is calculated and centered to obtain the centered frequency domain image. : ； In the formula, Represents the two-dimensional Fourier transform operator. This represents the Fourier transform centering operator; Let the frequency domain center column be , Represents the width of a two-dimensional astronomical image, preserving the column range. ,in The half-width parameter is used to determine the frequency range covered by the template and to construct the frequency domain mask. Make ； Frequency domain preservation is then performed to obtain the frequency domain image after frequency domain masking. : ； By performing inverse centering and inverse Fourier transform on it, a two-dimensional stripe template in the image domain can be obtained: ； in The operator for taking the real part is... This represents the two-dimensional inverse Fourier transform operator. This represents element-wise multiplication; These represent the row and column coordinates of the frequency domain image, respectively. Indicates the frequency domain mask in coordinates Pixel value at that location This represents the inverse centering operator of the Fourier transform.

9. The astronomical imaging fringe noise correction method based on low-frequency directional characteristics according to claim 1, characterized in that, In step S4, under the constraint of a fixed trapezoidal mask, the two-dimensional fringe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation and correction, and the image after fringe noise correction is output: ； in The template compensation strength coefficient is used to adjust the correction amplitude. To fix the trapezoidal mask, Represents a two-dimensional stripe template in the image domain. This represents a two-dimensional astronomical image after column-level brightness equalization.

10. An astronomical imaging stripe noise correction system based on low-frequency directional characteristics, characterized in that, The system includes: The vertical stripe brightness equalization module is used to perform vertical stripe brightness equalization on the acquired two-dimensional astronomical image data by using robust statistical estimation of the representative brightness of each column and performing column-level brightness shifting. This results in a column-level brightness equalized two-dimensional astronomical image. The low-frequency directional feature heatmap calculation module is used to slide sample on the two-dimensional astronomical image after column-level brightness equalization with a set window size and step size. For each window, it performs mean removal, windowing and two-dimensional Fourier transform to extract low-frequency ring energy and statistically analyze the directional energy distribution by angle binning. It calculates the low-frequency directional coefficient and weights the results of each window according to the coverage area to generate a low-frequency directional feature heatmap of the same size as the original image to characterize the spatial position and relative intensity of residual fringes. The trapezoidal fitting module is used to determine the contrast validity of the low-frequency directional feature heatmap and then perform threshold segmentation, morphological closing operation and denoising. Only the connected regions with the largest area and that meet the area ratio constraint are retained as the main region mask. The left and right boundary lines are fitted to the main region mask line by line, and the intersection of the two lines with the upper and lower boundaries of the image is used to obtain four vertices, generating a fixed trapezoidal mask. The stripe template estimation and compensation correction module is used to perform star point suppression or star point masking on the original two-dimensional astronomical image, and then perform a two-dimensional Fourier transform on the processed two-dimensional astronomical image to retain only the vertical frequency band near the center of the frequency domain to obtain the vertical stripes / low-frequency components. The inverse transform is used to obtain the two-dimensional stripe template in the image domain. Under the constraint of a fixed trapezoidal mask, the two-dimensional stripe template in the image domain is directly applied to the two-dimensional astronomical image after column-level brightness equalization for compensation correction, and the stripe noise corrected image is output.

Images