Method and system for detecting airglow gravity waves based on multi-stage filtering and deep learning

Abstract

The present application relates to the technical field of atmospheric optical remote sensing, and discloses a method and system for detecting airglow gravity waves based on multi-stage filtering and deep learning. The method comprises the following steps: after the star difference of adjacent airglow images is removed, the background trend is extracted by median filtering and Gaussian filtering in turn, and an enhanced image is obtained by geographic projection and adaptive stretching; after the enhanced image is added with a Kaiser window, a cover image is obtained by multi-scale Gabor filtering; a gravity wave frame is input into a convolutional neural network for detection; the current frame region is intercepted by an extended frame, and a windowed Fourier transform is performed to obtain an amplitude spectrum; meanwhile, the group velocity is calculated by using cross-correlation matching; the peak value with the smallest angle with the group velocity direction is selected in the amplitude spectrum to calculate the wavelength and direction; and the group velocity direction is corrected according to the relationship between the phase velocity direction and the group velocity direction. The present application solves the problem of automatic recognition of weak airglow signals, and realizes high-precision automatic detection of gravity waves and multi-parameter quantization.

Description

Technical Field

[0001] This invention relates to the field of atmospheric optical remote sensing technology, specifically to a method and system for detecting airglow gravity waves based on multi-level filtering and deep learning. Background Technology

[0002] Airglow imaging is an important tool for monitoring gravity wave activity in the middle and upper atmosphere. Typical airglow imagers acquire all-day nighttime images in the OH band with high temporal resolution. However, these images contain significant star interference, instrument-fixed mode noise, and a slowly changing static atmospheric background, making it extremely difficult to automatically identify and detect weak quasi-monochromatic gravity wave fluctuations. Furthermore, gravity wave trains are characterized by fine textures and low contrast, and these features are often difficult to directly identify due to imaging inhomogeneities and geographic projection distortion. Existing methods mostly rely on manual interpretation or simple frequency domain filtering, which is inefficient, highly subjective, and unsuitable for the batch processing of massive historical data.

[0003] In recent years, some studies have attempted to use machine learning methods for recognition and detection, but often the preprocessing is insufficient, resulting in excessive noise in the feature map, or the design is not tailored to the physical characteristics of gravity waves, and there is a lack of automated quantitative analysis of physical parameters such as wavelength and propagation direction, as well as intuitive multi-view comparison output. Summary of the Invention

[0004] The technical problem to be solved by this invention is: how to extract gasglow images (generally from a gasglow imager) from a gasglow imager. While removing strong interference such as star points and background in pixels, it enhances the gravity wave texture structure at different scales and uses deep learning to accurately identify the fluctuation region of gravity waves. At the same time, it automatically completes Fourier analysis and cross-correlation analysis to accurately invert the wavelength, propagation direction and group velocity of gravity wave activity, so as to meet the needs of automatic batch processing of large-scale observation data and scientific research.

[0005] To solve the above-mentioned technical problems, the present invention adopts the following technical solution: In a first aspect, the present invention provides a method for detecting airglow gravity waves based on multi-level filtering and deep learning, comprising: After performing median filtering on two adjacent geometrically corrected gasglow images, the difference between the two frames is calculated to obtain a difference image. The difference image is subjected to a second median filter and a Gaussian filter in sequence. The background trend extracted by the Gaussian filter is used for geographic projection, and then the enhanced image is obtained by adaptive contrast stretching. After applying a two-dimensional Kaiser window to the enhanced image, multi-scale and multi-directional Gabor filtering is performed. All filter responses are accumulated, normalized, and histogram equalized to obtain the Gabor coverage map. After copying the Gabor overlay map along the channel dimension into three channels, it is input into a pre-trained convolutional neural network model, which outputs a rectangular box of the gravity wave region. For each rectangle, expand it outward by a preset pixel width, extract the corresponding expanded region in the current frame enhanced image, window it, and perform a two-dimensional Fourier transform to obtain the amplitude spectrum. At the same time, use the current frame rectangle region as a template and the expanded region in the previous frame enhanced image as a search region to perform normalized cross-correlation matching to obtain the group velocity direction and magnitude. Select the peak value with the smallest angle between the phase velocity direction and the group velocity direction in the amplitude spectrum, and calculate the wavelength and propagation direction of the gravity wave within the rectangle. If the product of the phase velocity direction and the group velocity direction of the peak value is negative, then the group velocity direction is reversed.

[0006] In one embodiment, the standard deviation σ of the Gaussian filter is 2, which is used to specifically extract the background trend from the differential image after removing star points and impulse noise; the geographic projection is interpolated based on the background trend, and the adaptive contrast stretching includes: when the dynamic range of the geographically projected image is less than a set value, stretching is performed using the mean ±3σ method.

[0007] In one embodiment, the β parameter of the two-dimensional Kaiser window is 4, obtained by the outer product of two one-dimensional Kaiser windows, and is used to suppress the ringing effect generated by image edges during Gabor filtering; the multi-scale, multi-directional Gabor filtering includes 7 scales and 8 directions; the 8 directions are 0, π / 8, π / 4, 3π / 8, π / 2, 5π / 8, 3π / 4 and 7π / 8 respectively; the aspect ratio γ of the Gabor filter space is 0.5, the phase offset ψ is 0, and the standard deviation σ and wavelength λ of the filter are scaled proportionally with the kernel size.

[0008] In one embodiment, the convolutional neural network model employs the YOLOv8n architecture.

[0009] In one embodiment, the windowing uses a two-dimensional Hanning window; the amplitude spectrum is a logarithmic amplitude spectrum, and when searching for peaks in the amplitude spectrum, a central region of a set size is excluded to eliminate zero-frequency interference; the normalized cross-correlation is used to perform sub-pixel quadratic surface fitting on the matched response map to obtain the sub-pixel offset; The group velocity magnitude is calculated based on the subpixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the subpixel offset. When selecting the peak with the smallest angle, if the group velocity direction is valid, select the peak with the smallest angle between the phase velocity direction and the group velocity direction from the amplitude spectrum; otherwise, select the peak with the largest amplitude value.

[0010] In one embodiment, the group velocity magnitude is calculated based on the sub-pixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the sub-pixel offset, specifically including: ; ; For the group velocity magnitude, For the direction of group velocity, These are the sub-pixel offsets in the horizontal and vertical directions, respectively. For geographic resolution, This represents the time interval between consecutive frames.

[0011] Secondly, the present invention provides an airglow gravity wave detection system based on multi-level filtering and deep learning, comprising: The single-frame preprocessing module performs median filtering on two adjacent geometrically corrected gasglow images and then calculates the difference between the two frames to obtain a difference image. The inter-frame difference enhancement module performs a second median filter and a Gaussian filter on the difference image in sequence, uses the background trend extracted by the Gaussian filter for geographic projection, and then obtains the enhanced image through adaptive contrast stretching. The Gabor feature extraction module applies a two-dimensional Kaiser window to the enhanced image, performs multi-scale and multi-directional Gabor filtering, accumulates, normalizes, and histogram equalizes all filter responses to obtain a Gabor coverage map. The convolutional neural network detection module copies the Gabor cover map into three channels along the channel dimension, inputs it into a pre-trained convolutional neural network model, and outputs a rectangular box of the gravity wave region. The joint analysis module expands each rectangle outward by a preset pixel width, extracts the corresponding expanded region in the current frame's enhanced image, performs a two-dimensional Fourier transform after windowing to obtain the amplitude spectrum, and simultaneously uses the current frame's rectangle region as a template and the expanded region in the previous frame's enhanced image as a search region to perform normalized cross-correlation matching to obtain the group velocity direction and magnitude. The module then selects the peak with the smallest angle between the phase velocity direction and the group velocity direction in the amplitude spectrum, calculates the wavelength and propagation direction of the gravity wave within that rectangle, and reverses the group velocity direction if the product of the peak's phase velocity direction and the group velocity direction is negative.

[0012] In one embodiment, the standard deviation σ of the Gaussian filter is 2, which is used to specifically extract the background trend from the differential image after removing star points and impulse noise; the geographic projection is interpolated based on the background trend, and the adaptive contrast stretching includes: when the dynamic range of the geographically projected image is less than a set value, stretching is performed using the mean ±3σ method.

[0013] In one embodiment, the β parameter of the two-dimensional Kaiser window is 4, obtained by the outer product of two one-dimensional Kaiser windows, and is used to suppress the ringing effect generated by image edges during Gabor filtering; the multi-scale and multi-directional Gabor filtering includes 7 scales and 8 directions; the 8 directions are 0, π / 8, π / 4, 3π / 8, π / 2, 5π / 8, 3π / 4 and 7π / 8 respectively; the aspect ratio γ of the Gabor filter space is 0.5, the phase offset ψ is 0, and the standard deviation σ and wavelength λ of the filter are scaled proportionally with the kernel size; The convolutional neural network model uses the YOLOv8n architecture.

[0014] In one embodiment, the windowing uses a two-dimensional Hanning window; the amplitude spectrum is a logarithmic amplitude spectrum, and when searching for peaks in the amplitude spectrum, a central region of a set size is excluded to eliminate zero-frequency interference; the normalized cross-correlation is used to perform sub-pixel quadratic surface fitting on the matched response map to obtain the sub-pixel offset; The group velocity magnitude is calculated based on the subpixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the subpixel offset, specifically including: ; ; For the group velocity magnitude, For the direction of group velocity, These are the sub-pixel offsets in the horizontal and vertical directions, respectively. For geographic resolution, The time interval between consecutive frames; When selecting the peak with the smallest angle, if the group velocity direction is valid, select the peak with the smallest angle between the phase velocity direction and the group velocity direction from the amplitude spectrum; otherwise, select the peak with the largest amplitude value.

[0015] The system and method in this invention correspond to each other; the specific technical solutions applicable to the method are also applicable to the system.

[0016] Compared with the prior art, the beneficial technical effects of the present invention are: This invention combines Gaussian filtering for directional background trend extraction, Kaiser window and multi-scale Gabor filtering for noise reduction and enhancement, deep learning detection, and cross-correlation-spectrum joint dual-mode velocity analysis. It not only achieves fully automated high-precision detection of gravity wave regions, but also provides complete physical parameters such as phase velocity wavelength, direction, and group velocity magnitude and direction in a single process, providing an efficient and reliable solution for the batch scientific analysis of airglow historical data. Attached Figure Description

[0017] Figure 1This is the overall flowchart for automatic extraction of gravity wave structures and parameter calculation in OH gasglow images.

[0018] Figure 2 The images show a comparison of the effects of each step in the preprocessing and enhancement stages, where: (a) is the original geometrically corrected gasglow image; (b) is the image after median filtering for star removal and inter-frame differential enhancement; (c) is the noise residual image extracted during Gaussian smoothing; and (d) is the image after secondary median filtering and Gaussian smoothing.

[0019] Figure 3 The images show a comparison of the effects of the feature extraction stage, where: (a) is the enhanced image after geographic projection, adaptive contrast stretching, and the application of a Kaiser window; and (b) is the overlay map generated after multi-scale, multi-directional Gabor filtering.

[0020] Figure 4 The final output is a diagram showing the analysis results with latitude and longitude grids and gravity wave detection rectangles.

[0021] Figure 5 This is an F1 curve during the training process of the YOLOv8n model in an embodiment of the present invention.

[0022] Figure 6 This is a PR curve diagram during the training process of the YOLOv8n model in an embodiment of the present invention.

[0023] Figure 7 This is a flowchart of the method of the present invention. Detailed Implementation

[0024] A preferred embodiment of the present invention will now be described in detail with reference to the accompanying drawings.

[0025] like Figure 7 As shown, an airglow gravity wave detection method based on multi-level filtering and deep learning in this invention includes the following steps: S1. After performing median filtering to remove star points on two adjacent geometrically corrected gasglow images, calculate the difference between the two frames to obtain the difference image. S2, perform a second median filter and a Gaussian filter on the difference image in sequence, use the background trend extracted by the Gaussian filter for geographic projection, and then obtain the enhanced image through adaptive contrast stretching. S3. After applying a two-dimensional Kaiser window to the enhanced image, perform multi-scale and multi-directional Gabor filtering, accumulate, normalize and histogram equalize all filter responses to obtain the Gabor coverage map. S4, after copying the Gabor overlay map along the channel dimension into three channels, input it into a convolutional neural network model pre-trained with YOLOv8n, and output the rectangular box of the gravity wave region; S5. For each rectangle, expand it outward by a preset pixel width, extract the corresponding expanded region in the current frame enhanced image, window it, and perform a two-dimensional Fourier transform to obtain the amplitude spectrum. At the same time, use the current frame rectangle region as a template and the expanded region in the previous frame enhanced image as a search region to perform normalized cross-correlation matching to obtain the group velocity direction and magnitude. Select the peak value with the smallest angle between the phase velocity direction and the group velocity direction in the amplitude spectrum, calculate the wavelength and propagation direction of the gravity wave within the rectangle, and correct the group velocity direction based on the phase velocity direction of the peak value.

[0026] The present invention will be described in detail below in several parts.

[0027] 1. Single-frame preprocessing.

[0028] like Figure 1 As shown, the raw binary observation data stored in the OH gasglow imager is read, and negative pixels are replaced with preset background constants to remove bad pixels; the horizontal offset and rotation direction of the image are corrected to obtain a geometrically corrected gasglow image.

[0029] Specifically, the removal of bad pixels in single-frame preprocessing includes: replacing negative pixels in the original data with a preset saturation value backcon=20000; geometric correction includes pixel offset of shift1=9 and shift2=12 pixels according to the instrument calibration parameters (shift1, shift2, rotate1), rotating around the center by rotate1=34°, and filling the gaps with const1=1000.

[0030] 2. Inter-frame differential enhancement.

[0031] like Figure 1 As shown, large-size median filtering is applied to two adjacent frames of geometrically corrected airglow images, and overly bright star pixels are detected and replaced based on the global mean threshold to obtain a star-removed image. The difference between the current frame and the previous frame star-removed image is calculated to obtain a difference map to remove fixed noise. Small-size median filtering and Gaussian smoothing are applied to the difference map in sequence to obtain a smoothed difference map. Using a pre-established geographic projection coordinate system, the smoothed difference map is interpolated to an equidistant geographic grid to obtain a geographically corrected map. Numerical cropping and adaptive contrast stretching are performed on the geographically corrected map to generate an enhanced image.

[0032] Specifically, in the inter-frame differential enhancement, the large-size median filter window for star removal is 15×15 pixels, and the bright star threshold ratio is 10% of the global mean; in the secondary smoothing, the small-size median filter window is 5×5 pixels, and the Gaussian filter standard deviation σ=2.0. The purpose is not simply to smooth noise, but to specifically extract large-scale, slowly varying background trend components from the differential image after removing star points and impulse noise. In the geographic projection stage, Delaunay triangulation and a linear interpolator are used to resample the smoothed differential map from the instrument coordinates to a 1024×1024 grid in the range field1=[-343km,343km], and the pixel values ​​are truncated to [-800,800]. During adaptive contrast stretching, if the dynamic range is less than 100 pixel values, the mean ±3σ is used for cropping before linear stretching; otherwise, it is directly linearly stretched to 0-255.

[0033] 3. Gabor feature extraction.

[0034] After applying a two-dimensional Kaiser window to the enhanced image, multi-scale, multi-directional Gabor filtering is performed. The responses of each channel are summed, normalized, and histogram equalized to obtain the Gabor coverage map.

[0035] Specifically, when processing and enhancing the image, the image is first converted to floating-point type and the mean is subtracted, then multiplied by a two-dimensional Kaiser window with a β parameter value of 4. The attenuation of the Kaiser window at the image edges effectively prevents spectral leakage and ringing during Gabor convolution. Subsequently, a set of Gabor filters with 7 scales and 8 directions is applied. The kernel sizes for the 7 scales are scale1=31, scale2=41, scale3=51, scale4=61, scale5=71, scale6=81, and scale7=91, and the angles for the 8 directions are orient1=0, orient2=π / 8, orient3=π / 4, orient4=3π / 8, orient5=π / 2, orient6=5π / 8, orient7=3π / 4, and orient8=7π / 8. The filter parameters are γ=0.5, ψ=0, and σ and wavelength λ are scaled proportionally with the kernel size (when the reference σ=4.0, the reference λ=10.0, and the kernel size is 31). The 56 filter responses are summed, normalized to 0-255, and then histogram equalization is performed to obtain a high-contrast Gabor overlay map.

[0036] 4. Convolutional Neural Network Detection.

[0037] The Gabor overlay map is input into a pre-trained convolutional neural network (CNN) model, which outputs the coordinates of rectangular boxes representing suspected gravitational wave regions.

[0038] Specifically, the Gabor overlay image was input into a convolutional neural network model specifically trained for this type of data. The training dataset contained 656 images, using the YOLOv8n architecture, with an input image size of 640×640 pixels. Training lasted for 100 epochs with a batch size of 16. During training, the bounding box loss, classification loss, and distributed focal length (DFL) loss steadily decreased. On the validation set, the precision reached 0.9132, the recall was 0.8646, and the mAP50 was 0.9316. Because the input images had undergone fine-tuning and enhancement, the model did not need to learn complex interference patterns, significantly improving the detection rate and reducing false positives. The Gabor overlay image to be detected (single-channel grayscale) was copied into a three-channel image and input into the trained YOLOv8n model. The model output several rectangular boxes [x1, y1, x2, y2], each representing a suspected gravity wave region.

[0039] 5. Joint analysis of group velocity and phase velocity.

[0040] For each detection box, a preset number of pixels are expanded outwards. The corresponding expanded region in the current frame's enhanced image is extracted, windowed, and subjected to a two-dimensional Fast Fourier Transform (FFT) to obtain the logarithmic amplitude spectrum. A peak is searched within the amplitude spectrum (excluding the central region), but this peak exhibits 180° directional ambiguity. Simultaneously, using the original rectangular bounding box region of the current frame as a template, normalized cross-correlation matching is performed in the expanded region of the previous frame's enhanced image. Subpixel quadratic surface fitting is applied to the matching response map to obtain subpixel offsets. The magnitude and direction of the group velocity are calculated based on the offset, geographic resolution, and the time interval between consecutive frames. For candidate peaks in the amplitude spectrum, the peak with the smallest angle to the group velocity direction is selected. The wavelength is calculated from the spatial frequency of this peak, and its direction is the phase velocity propagation direction. Finally, the group velocity direction is corrected based on the sign relationship between the phase velocity direction and the group velocity direction.

[0041] In the joint analysis of group velocity and phase velocity, normalized cross-correlation matching adopts phase correlation or normalized cross-correlation methods and performs sub-pixel peak fitting. The inter-frame time difference and pixel displacement are combined with the geographic resolution per pixel to calculate the magnitude and direction of the group velocity. The preset pixel width for outward expansion is 100 pixels, and a two-dimensional Hanning window is used for windowing. The amplitude spectrum is a logarithmic amplitude spectrum, and the region with a center radius of 10 pixels is excluded when searching for peaks. When selecting the peak with the smallest angle with the group velocity direction, if there is no group velocity direction, the peak with the largest amplitude is directly used. The sign correction is: if the calculated phase velocity direction has the opposite sign to the group velocity direction, the group velocity direction is reversed.

[0042] Example: This embodiment uses int16 format time series images acquired by an OH gasglow imager with a field of view of 150° at a resolution of 1024×1024 pixels as the processing object. The original data consists of 940,009 .bin files, and after processing, 6,176 valid gasglow images are obtained.

[0043] Step 1: Single-frame preprocessing: The two-dimensional array `data` is read from the original binary file. Abnormally negative pixels caused by sensor saturation are replaced with a saturation value `backcon=20000` to avoid overexposure distortion. Pixels are then offset according to the instrument calibration parameters (`shift1`, `shift2`, `rotate1`) by 9 pixels (`shift1`) and 12 pixels (`shift2`), and rotated around the center by `rotate1` = 34°. This step eliminates bad pixels and geometric distortion, and unifies the image orientation.

[0044] Step 2, Inter-frame Differential Enhancement: Two adjacent frames of geometrically corrected airglow images are selected. First, a 15×15 pixel median filter is applied, with a threshold ratio set to 10% of the global mean. Only points where the brightness value changes beyond the threshold before and after filtering are replaced, resulting in a de-star image. The difference matrix is ​​obtained by subtracting the previous de-star image from the current one. A 5×5 pixel median filter and a Gaussian filter with σ=2.0 are then applied to the difference matrix to obtain a smoothed difference. Using a pre-generated Delaunay triangulation and linear interpolator, the smoothed difference is resampled from the instrument coordinates to a 1024×1024 grid within the range of [-343km, 343km], and the pixel values ​​are truncated to [-800, 800]. Then, adaptive contrast stretching is performed: if the dynamic range is less than 100 pixel values, the mean ± 3σ is used for cropping before linear stretching; otherwise, it is directly linearly stretched to 0-255.

[0045] Figure 2 This is a comparison chart showing the effects of each step in the preprocessing and enhancement stages. Figure 2 (a) in the image is the original geometrically corrected gasglow image; Figure 2 Image (b) is the image after median filtering for satellite removal and inter-frame difference enhancement; Figure 2 (c) in the figure is the noise residual map extracted during the Gaussian smoothing process; Figure 2 In the image, (d) represents the image after quadratic median filtering and Gaussian smoothing.

[0046] Step 3: Gabor Feature Extraction The enhanced image is converted to floating-point form and the mean is subtracted, then multiplied by a 2D Kaiser window with β=4 (obtained by the outer product of two 1D Kaiser windows). A set of Gabor filters with 7 scales (kernel sizes scale1=31, scale2=41, scale3=51, scale4=61, scale5=71, scale6=81, scale7=91) and 8 directions (orient1=0, orient2=π / 8, orient3=π / 4, orient4=3π / 8, orient5=π / 2, orient6=5π / 8, orient7=3π / 4, orient8=7π / 8) is then applied, with parameters γ=0.5, ψ=0, and σ and λ scaled proportionally with the kernel size (reference σ=4.0, reference λ=10.0, kernel size 31). The 56 filter responses are accumulated, normalized to 0-255, and then histogram equalization is performed to obtain a high-contrast Gabor overlay map.

[0047] Figure 3 This is a comparison chart showing the results of the feature extraction stage. Figure 3 (a) in the image is the enhanced image after geographic projection, adaptive contrast stretching, and the application of a Kaiser window; Figure 3 (b) in the image is the Gabor coverage map generated after multi-scale, multi-directional Gabor filtering.

[0048] Step 4: Convolutional Neural Network Detection: First, the Gabor overlay images were manually labeled to select gravitational wave regions, constructing a dataset of 656 images. The YOLOv8n architecture was used, with input images uniformly scaled to 640×640 pixels. Training was performed for 100 epochs with a batch size of 16 and an initial learning rate of 0.001. After training, the precision was 0.9132, recall 0.8646, and mean precision (mAP50) 0.9316 on the validation set. Since the convolutional neural network model has three input channels, the single-channel Gabor overlay grayscale image was copied three times along the channel dimension and merged into a three-channel image, ensuring identical pixel values ​​in each channel. This three-channel image was then input into the pre-trained convolutional neural network model for gravitational wave region detection, outputting several bounding boxes [x1, y1, x2, y2] and confidence scores.

[0049] Figure 5 This is the F1 curve (BoxF1_curve) during the training process of the YOLOv8n model in this embodiment of the invention. The horizontal axis is the confidence threshold and the vertical axis is the F1 score, which is used to evaluate the overall performance of the model under different thresholds.

[0050] Figure 6 This is a BoxPR curve (BoxPR_curve) during the training process of the YOLOv8n model in this embodiment of the invention. The horizontal axis represents recall, the vertical axis represents precision, the area under the curve is the mean precision (AP), and the average of each class is mAP.

[0051] Step 5: Joint analysis of group velocity and phase velocity: For each detection box, expand outward by 100 pixels (without exceeding the image boundary). Extract the expanded region from the current frame's enhanced image, multiply by a two-dimensional Hanning window, and perform a two-dimensional Fast Fourier Transform (FFT) to obtain the logarithmic amplitude spectrum log(1+|F|). After masking the center radius by 10 pixels (to eliminate zero-frequency interference), find the spatial frequency corresponding to the local peak. The direction corresponding to this peak has a 180° ambiguity. Simultaneously, using the original rectangular bounding box region of the current frame as a template, normalized cross-correlation matching is performed in the extended region of the enhanced image from the previous frame. Subpixel quadratic surface fitting is then applied to the matching response map to obtain the subpixel offset. The magnitude of the group velocity is calculated based on the offset, geographic resolution r = 0.67 km / pixel, and the time interval Δt between consecutive frames. and direction The unit of group velocity is For each candidate peak in the amplitude spectrum, calculate the spatial orientation angle corresponding to each peak. Choose the direction of the group velocity The candidate peak with the smallest included angle is taken as the phase velocity direction. If there is no group velocity direction, then select the peak value with the largest amplitude. The spatial frequency of this peak value is then used. Wavelength obtained by conversion (km), its direction This refers to the direction of phase velocity propagation (in degrees, measured counterclockwise from due east). Finally, if the phase velocity direction... Group velocity direction If the signs are opposite, then the group velocity direction will be reversed. Invert.

[0052] Step Six: Batch Processing and Result Output: Figure 4 Display the Gabor overlay map or the original enhanced image, overlaid with rectangular boxes and numbers 1, 2, and 3; Table 1 shows... Figure 4 The parameter table attached to the analysis results graph, Table 1, lists the wavelength (in km), phase velocity direction (in degrees), group velocity (in m / s), and group velocity direction (in degrees) for each box. All steps can be executed in parallel on multiple date data, enabling unattended batch analysis.

[0053] Table 1

[0054] The terminology used herein is for the purpose of describing particular embodiments only and is not intended to limit the invention. The terms “comprising,” “including,” etc., as used herein indicate the presence of the stated features, steps, operations, and / or components, but do not exclude the presence or addition of one or more other features, steps, operations, or components.

[0055] It should be understood that although the steps in the flowcharts of the accompanying drawings are shown sequentially as indicated by the arrows, these steps are not necessarily executed in the order indicated by the arrows. Unless explicitly stated herein, there is no strict order restriction on the execution of these steps, and they can be executed in other orders. Moreover, at least some of the steps in the flowcharts of the accompanying drawings may include multiple steps or stages, which are not necessarily completed at the same time, but may be executed at different times, and the execution order of these steps or stages is not necessarily sequential, but may be performed alternately or in turn with other steps or at least some of the steps or stages of other steps.

[0056] Based on the description of the above method embodiments, the present invention also provides a system. The system may be a system that uses software (applications), modules, components, servers, clients, etc., using the methods described in the embodiments of this specification, combined with necessary implementation hardware. Since the implementation schemes and methods for solving the problem are similar, the specific system implementations in the embodiments of this specification can be found in the implementations of the foregoing methods, and repeated details will not be described again. Although the system is preferably implemented in software, hardware implementation, or a combination of software and hardware, is also possible and contemplated.

[0057] The technical features of the above embodiments can be combined in any way. For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.

[0058] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the invention can be implemented in other specific forms without departing from its spirit or essential characteristics. Therefore, the embodiments should be considered in all respects as exemplary and non-limiting, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of equivalents of the claims are intended to be included within the present invention, and no reference numerals in the claims should be construed as limiting the scope of the claims.

[0059] Furthermore, it should be understood that although this specification describes embodiments, not every embodiment contains only one independent technical solution. This narrative style is merely for clarity. Those skilled in the art should consider the specification as a whole, and the technical solutions in each embodiment can also be appropriately combined to form other embodiments that can be understood by those skilled in the art.

Claims

1. A method for detecting airglow gravity waves based on multi-level filtering and deep learning, characterized in that, include: After performing median filtering on two adjacent geometrically corrected gasglow images, the difference between the two frames is calculated to obtain a difference image. The difference image is subjected to a second median filter and a Gaussian filter in sequence. The background trend extracted by the Gaussian filter is used for geographic projection, and then the enhanced image is obtained by adaptive contrast stretching. After applying a two-dimensional Kaiser window to the enhanced image, multi-scale and multi-directional Gabor filtering is performed. All filter responses are accumulated, normalized, and histogram equalized to obtain the Gabor coverage map. After copying the Gabor overlay map along the channel dimension into three channels, it is input into a pre-trained convolutional neural network model, which outputs a rectangular box of the gravity wave region. For each rectangle, expand it outward by a preset pixel width, extract the corresponding expanded region in the current frame enhanced image, window it, and perform a two-dimensional Fourier transform to obtain the amplitude spectrum. At the same time, use the current frame rectangle region as a template and the expanded region in the previous frame enhanced image as a search region to perform normalized cross-correlation matching to obtain the group velocity direction and magnitude. Select the peak with the smallest angle between the phase velocity direction and the group velocity direction in the amplitude spectrum, calculate the wavelength and propagation direction of the gravity wave within the rectangle, and if the product of the phase velocity direction and the group velocity direction of the peak is negative, then reverse the group velocity direction. The windowing uses a two-dimensional Hanning window; the amplitude spectrum is a logarithmic amplitude spectrum, and when searching for peaks in the amplitude spectrum, a central region of a set size is excluded to eliminate zero-frequency interference; the normalized cross-correlation is used to perform sub-pixel quadratic surface fitting on the matched response map to obtain the sub-pixel offset; The group velocity magnitude is calculated based on the subpixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the subpixel offset. When selecting the peak with the smallest angle, if the group velocity direction is valid, select the peak with the smallest angle between the phase velocity direction and the group velocity direction from the amplitude spectrum; otherwise, select the peak with the largest amplitude value.

2. The airglow gravity wave detection method based on multi-level filtering and deep learning according to claim 1, characterized in that, The standard deviation σ=2 of the Gaussian filter is used to specifically extract background trends from the difference image after removing star points and impulse noise; The geographic projection is interpolated based on the background trend, and the adaptive contrast stretching includes: when the dynamic range of the geographically projected image is less than a set value, stretching is performed using the mean ±3σ method.

3. The airglow gravity wave detection method based on multi-level filtering and deep learning according to claim 1, characterized in that, The β parameter of the two-dimensional Kaiser window is 4, obtained by the outer product of two one-dimensional Kaiser windows, and is used to suppress the ringing effect generated by image edges during Gabor filtering. The multi-scale and multi-directional Gabor filtering includes 7 scales and 8 directions. The 8 directions are 0, π / 8, π / 4, 3π / 8, π / 2, 5π / 8, 3π / 4 and 7π / 8, respectively. The aspect ratio of the Gabor filter space is γ=0.5, the phase offset is ψ=0, and the standard deviation σ and wavelength λ of the filter are scaled proportionally with the kernel size.

4. The airglow gravity wave detection method based on multi-level filtering and deep learning according to claim 1, characterized in that, The convolutional neural network model uses the YOLOv8n architecture.

5. The airglow gravity wave detection method based on multi-level filtering and deep learning according to claim 1, characterized in that, The group velocity magnitude is calculated based on the sub-pixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the sub-pixel offset, specifically including: ； ； For the group velocity magnitude, For the direction of group velocity, These are the sub-pixel offsets in the horizontal and vertical directions, respectively. For geographic resolution, This represents the time interval between consecutive frames.

6. A system for detecting airglow gravity waves based on multi-level filtering and deep learning, characterized in that, include: The single-frame preprocessing module performs median filtering on two adjacent geometrically corrected gasglow images and then calculates the difference between the two frames to obtain a difference image. The inter-frame difference enhancement module performs a second median filter and a Gaussian filter on the difference image in sequence, uses the background trend extracted by the Gaussian filter for geographic projection, and then obtains the enhanced image through adaptive contrast stretching. The Gabor feature extraction module applies a two-dimensional Kaiser window to the enhanced image, performs multi-scale and multi-directional Gabor filtering, accumulates, normalizes, and histogram equalizes all filter responses to obtain a Gabor coverage map. The convolutional neural network detection module copies the Gabor cover map into three channels along the channel dimension, inputs it into a pre-trained convolutional neural network model, and outputs a rectangular box of the gravity wave region. The joint analysis module expands each rectangle outward by a preset pixel width, extracts the corresponding expanded region in the current frame's enhanced image, performs a two-dimensional Fourier transform after windowing to obtain the amplitude spectrum, and simultaneously uses the current frame's rectangle region as a template and the expanded region in the previous frame's enhanced image as a search region to perform normalized cross-correlation matching to obtain the group velocity direction and magnitude; selects the peak with the smallest angle between the phase velocity direction and the group velocity direction in the amplitude spectrum, calculates the wavelength and propagation direction of the gravity wave within that rectangle, and if the product of the peak's phase velocity direction and the group velocity direction is negative, then the group velocity direction is reversed; The windowing uses a two-dimensional Hanning window; the amplitude spectrum is a logarithmic amplitude spectrum, and when searching for peaks in the amplitude spectrum, a central region of a set size is excluded to eliminate zero-frequency interference; the normalized cross-correlation is used to perform sub-pixel quadratic surface fitting on the matched response map to obtain the sub-pixel offset; The group velocity magnitude is calculated based on the subpixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the subpixel offset. When selecting the peak with the smallest angle, if the group velocity direction is valid, select the peak with the smallest angle between the phase velocity direction and the group velocity direction from the amplitude spectrum; otherwise, select the peak with the largest amplitude value.

7. The airglow gravity wave detection system based on multi-level filtering and deep learning according to claim 6, characterized in that, The standard deviation σ=2 of the Gaussian filter is used to specifically extract background trends from the difference image after removing star points and impulse noise; The geographic projection is interpolated based on the background trend, and the adaptive contrast stretching includes: when the dynamic range of the geographically projected image is less than a set value, stretching is performed using the mean ±3σ method.

8. The airglow gravity wave detection system based on multi-level filtering and deep learning according to claim 6, characterized in that, The β parameter of the two-dimensional Kaiser window is 4, obtained by the outer product of two one-dimensional Kaiser windows, and is used to suppress the ringing effect generated by image edges during Gabor filtering; the multi-scale and multi-directional Gabor filtering includes 7 scales and 8 directions; the 8 directions are 0, π / 8, π / 4, 3π / 8, π / 2, 5π / 8, 3π / 4 and 7π / 8 respectively; the aspect ratio γ of the Gabor filter space is 0.5, the phase offset ψ is 0, and the standard deviation σ and wavelength λ of the filter are scaled proportionally with the kernel size; The convolutional neural network model uses the YOLOv8n architecture.

9. The airglow gravity wave detection system based on multi-level filtering and deep learning according to claim 6, characterized in that, The group velocity magnitude is calculated based on the sub-pixel offset, geographic resolution, and the time interval between consecutive frames. The group velocity direction is the arctangent value in the opposite direction of the sub-pixel offset, specifically including: ； ； For the group velocity magnitude, For the direction of group velocity, These are the sub-pixel offsets in the horizontal and vertical directions, respectively. For geographic resolution, This represents the time interval between consecutive frames.

Images