A runway pavement flatness identification method based on flight data

Abstract

The present application belongs to the runway pavement flatness evaluation technical field, specifically relates to a runway pavement flatness identification method based on flight data, and aims to solve the problems that traditional runway pavement detection is difficult to popularize and difficult to meet the dynamic monitoring demand of pavement state for aviation safety. The present application is aimed at the trend change of QAR characteristic parameters of the airplane in the take-off or landing taxiing stage, uses a time series anomaly detection method to judge whether there is an abnormal change trend, so as to reflect whether there is a pavement uneven condition; the acceleration auxiliary interpolation method is used to calculate the uneven interval, and the position of the uneven zone is positioned by using the longitude and latitude of the airborne GPS, and the distance from the runway entrance. The present application changes the airplane itself into a mobile detection sensor, upgrades from periodic static mapping to dynamic monitoring along with the flight, and finally realizes the acquisition of pavement real running state data with low cost and high timeliness.

Description

Technical Field

[0001] This invention belongs to the field of runway pavement smoothness assessment technology, and specifically relates to a runway pavement smoothness identification method based on flight data. Background Technology

[0002] Vertical load factor (VRTG), a core parameter measuring the ratio of aircraft lift to gravity, is primarily manifested as ground support force during the aircraft's ground taxiing phase. Under static or ideal planar motion conditions, the ground support force and gravity are balanced, and the VRTG value stabilizes at 1g. However, when an aircraft taxis over uneven or sloping runways, the ground support force dynamically changes with speed, trajectory radius, and runway surface undulations, leading to significant fluctuations in VRTG. This parameter variation is closely related to runway surface conditions, providing crucial data for runway surface inspection.

[0003] The Quick Access Recorder (QAR), as a key airborne flight data recorder, collects and stores over 2,000 flight parameters in real time through the Flight Data Acquisition Unit (FDAU). These parameters cover multiple dimensions of aircraft operation, and its ability to quickly extract raw data provides an efficient and rich data foundation for subsequent in-depth analysis, making it an important data source for aircraft operational status research.

[0004] The takeoff and landing roll phases are key aspects of aviation safety. The takeoff roll begins with the aircraft entering and aligning with the runway, engaging takeoff thrust, and ending with the main wheels leaving the ground. The landing roll begins with the main wheels touching down and ends with the tail section leaving the runway. During these two phases, the aircraft is in close contact with the runway, and the runway surface condition directly affects the aircraft's operational safety and performance; therefore, accurate detection of the runway surface condition is crucial.

[0005] Due to varying runway foundation conditions, and the possibility of runway surface subsidence during airport operations, the actual aircraft operating experience may differ from the information published in airport regulations. From the user's perspective, airlines cannot directly and quickly obtain a visual understanding of the runway surface. When aircraft encounter uneven surfaces during high-speed taxiing, it increases the fatigue risk to the landing gear structure, interferes with pilot observation and judgment, and affects the passenger experience.

[0006] Currently, traditional methods for detecting runway pavement slope and smoothness, such as using professional instruments like laser rangefinders or satellite remote sensing technology, have significant limitations: these methods are not only expensive to purchase and maintain, and have complex operating procedures, making them difficult to implement in small and medium-sized airports; more importantly, they can only obtain static data on the geometric shape of the pavement, which is insufficient to meet the dynamic monitoring needs of aviation safety for pavement condition. Summary of the Invention

[0007] To address the aforementioned problems in existing technologies, namely that traditional runway pavement inspection relies on specialized instruments or satellite remote sensing technology and can only acquire static geometric data of the pavement, resulting in its limited widespread application and inability to meet the dynamic monitoring requirements of aviation safety for pavement condition, this invention, in its first aspect, proposes a runway pavement smoothness identification method based on flight data. This method includes the following steps:

[0008] S1, obtain the QAR dataset of the target aircraft during the takeoff and landing taxiing phase of the first runway of the target airport and extract the VRTG parameters to generate a time series set X1 of VRTG parameters; the parameters extracted from the QAR dataset also include ground speed, longitudinal overload, and vertical overload; the first runway is any runway of the target airport; S2, perform first-order difference processing on X1 to obtain the difference time series X2; S3, perform anomaly detection on X2 to obtain abnormal data of VRTG parameters; Based on the waveform changes of the abnormal data of the VRTG parameters, the vibration pattern is obtained, and then the uneven area of ​​the first runway surface is obtained. S4. Calculate the horizontal acceleration using the longitudinal overload of the target aircraft, and then use acceleration-assisted interpolation to complete the frequency of the ground speed parameter of the target aircraft to the frequency of the VRTG parameter. Based on the ground speed parameters after frequency completion, the length of the pit in the uneven section of the first runway surface is obtained by ground speed integration calculation. The vertical velocity is obtained by processing the vertical overload, and then the equivalent depth of the pit in the uneven section of the first runway surface is calculated according to the energy conservation model. S5, based on the length, locate the position of the pit within the uneven section of the first runway surface.

[0009] The beneficial effects of this invention are: 1) By analyzing the trend changes of QAR characteristic parameters during takeoff or landing rollouts, a time-series anomaly detection method is used to determine whether there are abnormal trends, thereby indicating the presence of pavement irregularities. Acceleration-assisted interpolation is used to calculate the irregularity intervals, and the location of the irregular area and its distance from the runway threshold is determined using the onboard GPS latitude and longitude recorded by the aircraft's QAR. By turning the aircraft itself into a mobile detection sensor, it is possible to determine whether the pavement is smooth and the distribution of irregularities. This upgrades from periodic static mapping to dynamic monitoring that occurs with flight operations, ultimately achieving a technological breakthrough in obtaining real pavement operational status data at low cost and high timeliness, providing a new dimension for risk prevention and control in the field of aviation safety. 2) Unlike instrument measurements, QAR data provides a direct reflection of the impact of pavement irregularities on aircraft, enabling efficient identification, timely assessment, and the implementation of appropriate measures to reduce operational risks through the identified flight status. Attached Figure Description

[0010] Other features, objects, and advantages of this application will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings: Figure 1 This is a flowchart of the steps of a runway pavement smoothness identification method based on flight data according to the present invention.

[0011] Figure 2 This invention relates to a method for identifying runway pavement smoothness based on flight data, specifically the VRTG parameter waveform diagram of runway 07 at Ganzhou Huangjin Airport during the takeoff and taxiing phase.

[0012] Figure 3 This invention relates to a method for identifying runway pavement smoothness based on flight data, specifically the VRTG parameters and related parameter waveforms for the takeoff and taxiing phase of runway 07 at Ganzhou Huangjin Airport.

[0013] Figure 4 This invention relates to a method for identifying runway pavement smoothness based on flight data, specifically a waveform diagram of anomalies in VRTG parameters during the takeoff and taxiing phase of runway 07 at Ganzhou Huangjin Airport.

[0014] Figure 5 This invention relates to a method for identifying runway pavement smoothness based on flight data, specifically the vibration waveform diagram of VRTG parameters during the takeoff and taxiing phase of runway 07 at Ganzhou Huangjin Airport.

[0015] Figure 6 This is a waveform diagram showing the correspondence between low-speed parameters and VRTG parameters in a runway pavement smoothness identification method based on flight data according to the present invention.

[0016] Figure 7 This is a coordinate point diagram of the maximum deviation of VRTG parameters during the takeoff and taxiing phase of runway 07 at Ganzhou Huangjin Airport, based on the runway pavement smoothness identification method of the present invention.

[0017] Figure 8 This invention relates to an electronic map of Ganzhou Huangjin Airport, which is based on a method for identifying runway pavement smoothness using flight data. Detailed Implementation

[0018] The present application will now be described in further detail with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are for illustrative purposes only and are not intended to limit the invention. Furthermore, it should be noted that, for ease of description, only the parts relevant to the invention are shown in the accompanying drawings.

[0019] It should be noted that, unless otherwise specified, the embodiments and features described in this application can be combined with each other. This application will now be described in detail with reference to the accompanying drawings and embodiments.

[0020] To more clearly illustrate the runway pavement smoothness identification method based on flight data of the present invention, the following will be combined with... Figures 1 to 8 The steps in the embodiments of the present invention will be described in detail below.

[0021] The first embodiment of the present invention provides a method for runway pavement smoothness identification based on flight data. (See also...) Figure 1 This includes the following steps: S1, acquire the QAR dataset of the target aircraft during the takeoff and landing taxiing phase (i.e., the takeoff and landing taxiing phases) of the first runway of the target airport and extract the VRTG parameters (i.e., vertical overload parameters), and generate a time series set X1 of VRTG parameters; the parameters extracted from the QAR dataset also include ground speed, longitudinal overload, and vertical overload; the first runway is any runway of the target airport. In this embodiment, the QAR dataset of N segments during the takeoff and landing taxiing phase of the first runway of the target airport is obtained and VRTG parameters are extracted, where N≥1; If N > 1, the extracted VRTG parameters are superimposed. Example of multi-segment overlay analysis: Data from 10 flight segments departing from runway 07 at Ganzhou Huangjin Airport were selected; Extract VRTG parameters for the takeoff roll phase of each flight segment (from the engagement of takeoff thrust to the main wheels leaving the ground); The VRTG time series of 10 flight segments were superimposed and averaged to obtain the comprehensive fluctuation curve; When the aircraft reached a certain point during its runway roll, significant fluctuations in VRTG (range 0.6-1.5g) occurred, reflecting a dynamic imbalance between ground support force and gravity, indicating an unevenness at the corresponding location on the runway (see...). Figure 2 ); Example of multi-parameter correlation analysis for a single flight segment: Select B-19 Data for one flight segment from runway 07 at Ganzhou Huangjin Airport; Extract VRTG parameters during the takeoff roll phase (from the engagement of takeoff thrust to the main wheels leaving the ground) of this flight segment. At the same time, extract relevant parameters, including the aircraft air-to-ground switch, throttle stick angle, airspeed, ground speed, pitch angle, heading, captain / first officer's bilateral pitch control, vertical overload, rate of climb, and radio altitude. Observe the changes in other parameters when VRTG fluctuates. When the airspeed reaches 100 knots (UTC 12:53:35), the VRTG begins to fluctuate (0.75-1.42g) for about 2 seconds. During the fluctuation, the throttle lever angle and pitch control remained unchanged, ruling out human intervention. At 12:53:44 UTC, the crew pulled back the control stick and the aircraft took off normally, confirming that the fluctuations were caused by pavement factors (see...). Figure 3 ); S2, perform first-order difference processing on X1 to obtain the difference time series X2; In this embodiment, X1 is subjected to first-order differencing to obtain the difference time series X2. The method is as follows: Calculate the first difference of each element in the time series set X1 to obtain the differenced time series X2; the formula is as follows: In the formula, It is a difference operator; for Each element in ; The difference sequence reflects the changes in VRTG parameters and is used to capture the abrupt changes caused by pavement irregularities; The following is a physical explanation for why a first-order difference is necessary: ​​based on the impulse theorem... F Δ t = m Δ v When the plane is at high speed v The gliding distance is L When the pavement is in an uneven area, the contact time between the wheel and the pavement is: Δ t = vL The impact force generated at this time F With Δ t Inversely proportional to speed v Proportional, when Δ t In a very short time (such as when there are sharp protrusions or depressions on the pavement), the impact force increases significantly, which is manifested as a sharp rise in the VRTG parameter. At the same time, because the damping effect of the aircraft's shock absorption system cannot absorb the energy completely in an instant, the wheels rebound, causing the VRTG to drop rapidly in a short period of time, forming a characteristic "peak-valley" fluctuation pattern. S3, perform anomaly detection on X2 to obtain abnormal data of VRTG parameters; In this embodiment, the abnormal data of the VRTG parameter is obtained through X2 anomaly detection, and the method is as follows: Calculate the probability density P of each difference value in X2, and fit the probability density statistical distribution using the kernel density estimation method (KDE); Each probability density P is compared with a preset probability density threshold P0. If P≤P0, the VRTG parameter of the next time point in X1 corresponding to the difference value is determined to be abnormal data, thereby obtaining an abnormal data set. The preset probability density threshold P0 is determined by selecting a skewness range of 3% based on the distribution characteristics of historical data. The calculated value of the threshold P0 includes -0.926. Based on the waveform changes of the abnormal data of the VRTG parameters, the vibration pattern is obtained, and then the uneven area of ​​the first runway surface is obtained. In this embodiment, the deviation degree of each abnormal point is calculated to extract the point corresponding to the maximum deviation degree as a quantitative indicator of the severity of pavement unevenness. The deviation degree is calculated as follows: ,in: This represents the probability density value of a single data point. This represents the mean of the probability density of normal samples. The standard deviation represents the probability density of normal samples; The outlier dataset includes a sequence number, index, VRTG difference value, probability density P, and deviation degree D for each outlier, as shown in Table 1. Table 1 lists the sequence number, index, VRTG difference value, probability density P, and deviation degree D for 13 outliers. Table 1:

[0022] If outliers appear consecutively in the outlier dataset and exceed a set duration, the segment is determined to be an abnormal section caused by pavement unevenness. See [link to relevant documentation]. Figure 4 ; The duration includes the abnormal duration calculated based on the sampling frequency of the VRTG parameters and the maximum index difference of the outlier. Based on the index numbers in the table above, determine the continuity and duration of the anomalies. Since the VRTG parameter has a frequency of 8Hz, the difference between the maximum and minimum index numbers of the 13 anomalies is 20. Therefore, the anomalies occurred within 3 seconds, indicating that the anomalies were continuous. The following analysis of the VRTG waveform vibration pattern illustrates the vibration law of VRTG on uneven pavement surfaces. When an aircraft wheel passes over a pavement pothole, the shock absorber strut undergoes a dynamic deformation process: Entering the pothole stage: Ground support force decreases sharply, the shock absorber strut extends, and the VRTG value decreases; At the bottom of the pothole stage: After the wheel touches the bottom, the shock absorber strut compresses rapidly, and the VRTG value rises rapidly; Leaving the pothole stage: The shock absorber strut rebounds due to the release of elastic energy, and the VRTG value decreases again, subsequently decaying to a stable state according to the damped vibration law. This process results in VRTG exhibiting a characteristic "decreasing-increasing-decreasing" change pattern (see...). Figure 5In a real-world test case at runway 07 of Ganzhou Airport: VRTG vibration characteristics during takeoff taxiing: Vibration frequency: approximately 2 Hz (period 0.5~0.7 seconds); Initial amplitude: approximately 0.2g; Normal amplitude: According to the regular oscillation period trend, VRTG at the blue arrow should be the peak, but VRTG shows a significant decrease; Abnormal amplitude: It expands to about 0.7g at the red dotted line and lasts for 3 cycles; Special phenomenon: When the aircraft's takeoff speed approaches VR (nose wheel lift-off speed), unevenness of the runway may cause the nose wheel to lift off the ground briefly (as shown by the second blue arrow), which manifests as increased VRTG fluctuations; Based on the vibration patterns of the VRTG waveform, the starting point of the anomaly was determined, thus identifying the uneven section of the first runway surface. According to the VRTG vibration analysis, the following was determined: Figure 5 The abnormal range is 11-17 (corresponding to 0.875 seconds); The first half (11-14): The amplitude changes relatively gently, corresponding to the gradual change area at the entrance of the pit (28.7 meters in length). The latter half (15-17): The amplitude increases significantly, corresponding to the area of ​​abrupt change in pit depth (length 14.5 meters). Based on the waveform changes of the abnormal data of the VRTG parameters, the vibration pattern is obtained, and then the unevenness range of the first runway surface is determined. Further improvements are made, and the unevenness range includes three types: long-wave undulation unevenness, mixed unevenness, and short-wave impact unevenness. The method is as follows: S3.1 Calculate the degree of deviation for each abnormal data in the abnormal data set (calculation method as above); the degree of deviation is the standardized deviation of the probability density value corresponding to the abnormal data relative to the statistical distribution of the probability density of normal samples; S3.2, Perform multi-resolution wavelet decomposition on the original VRTG parameter time series corresponding to the abnormal data set to extract low-frequency feature components and high-frequency feature components of a set Hertz (e.g., low-frequency feature components of 1-4Hz and high-frequency feature components of 8-16Hz). S3.3, Calculate the proportion of the energy of the low-frequency characteristic component in the total characteristic energy. , ,in, , These represent the signal energy of the corresponding frequency band. When the proportion is greater than the first threshold (preferably 0.7), it is determined to be a long-wave undulation type unevenness; when the proportion is between the first threshold and the second threshold (preferably 0.3), it is determined to be a mixed type unevenness; when the proportion is less than the second threshold, it is determined to be a short-wave impulse type unevenness. The initial time length is determined according to the determined unevenness type (specifically, a, the initial spatial length Linit is determined according to the determined unevenness type: long-wave undulation type unevenness: Linit = 3.0m; mixed type unevenness: Linit = 1.5m; short-wave impulse type unevenness: Linit = 0.8m; b, the initial spatial length Linit is combined with the current ground speed v to calculate the initial time length Tinit, Tinit = Linit / v). S3.4, taking the abnormal data with the largest deviation as the center, extend to both sides by half of the initial time length, i.e., Tinit / 2, to form the initial abnormal interval; S3.5 Calculate the comprehensive anomaly index for each time point within the initial anomaly interval. The comprehensive anomaly index is a weighted combination value based on the degree of deviation, the rate of change of the VRTG parameter, and local volatility. The weighting coefficients of the weighted combination are adaptively adjusted according to the unevenness type. Specifically, the expression for calculating the comprehensive anomaly index is: ; in, The degree of deviation calculated in step S3.1 All are weighting coefficients. In this embodiment, the coefficients are (0.4, 0.3, 0.3) for the long-wave undulation type, (0.5, 0.3, 0.2) for the mixed type, and (0.6, 0.2, 0.2) for the short-wave impulse type. It is the absolute value (i.e., the rate of change) of the first difference of the VRTG parameter. The local standard deviation is calculated as follows: the standard deviation of the VRTG parameter centered at the current sampling point t, with two sampling points before and after it (a total of 5 points); S3.6, Starting from the boundary of the initial abnormal interval, iteratively expand outwards, and determine whether to include adjacent time points in the uneven interval based on the comparison between the comprehensive abnormal index and a preset threshold. When the comprehensive abnormal index of multiple consecutive time points (e.g., 3) is lower than the preset threshold... (Preferred) Expansion stops when the value reaches 0.5, the preset threshold. The index is dynamically determined based on the statistical characteristics of the comprehensive anomaly index within the initial anomaly interval. ,in, This represents the average value of the comprehensive anomaly index within the initial anomaly interval. , The population standard deviation of the comprehensive anomaly index within the initial anomaly interval. Ω is the set of time points within the initial abnormal interval, and N is the total number of time points within the interval. S3.7, the time boundary after iterative expansion is determined as the final uneven interval. and output the corresponding unevenness type; Spatial domain verification: Total length of the crater calculated by the 8Hz ground velocity integral: 43.2 meters; The ground speed estimation result at 1Hz is 50 meters (error 16%). This difference verifies the effectiveness of high-frequency interpolation in improving positioning accuracy. Based on this, when calculating the length and depth of the pits in the uneven section of the first runway surface, the ground speed parameters and vertical velocity are supplemented by frequency. S4. Calculate the horizontal acceleration using the longitudinal overload of the target aircraft, and then use acceleration-assisted interpolation to complete the frequency of the ground speed parameter of the target aircraft to the frequency of the VRTG parameter. In this embodiment, the frequency of the ground speed parameter of the target aircraft is supplemented to the frequency of the VRTG parameter using an acceleration-assisted interpolation method. The method is as follows: Based on the initial velocity and current horizontal acceleration of the ground speed parameters, interpolation calculations are performed according to a preset time interval. The sampling interval of the target aircraft is divided into multiple interpolation points for interpolation completion, so that the sampling frequency of the ground speed parameters reaches the frequency of the VRTG parameters, thus obtaining the frequency-completed ground speed parameters. That is, the high-frequency ground velocity sequence is obtained, and the interpolation formula expression is: ; in: For the first The known ground velocity value at each moment; For the first The horizontal acceleration calculated from the longitudinal overload at each moment , The acceleration due to gravity is 9.80665 m / s²; Δt represents the time interval, and the sampling interval N of the target aircraft is divided into N / M interpolation points for interpolation completion, where M is the frequency of the VRTG parameter; Based on the ground speed parameters after frequency completion, the length of the pit in the uneven section of the first runway surface is obtained by ground speed integration calculation, and the expression for the length L is: ; The following detailed explanation, with examples, illustrates the acquisition of ground speed and longitudinal acceleration parameters. When the aircraft pitch angle is 0, the longitudinal acceleration is determined to be equal to the horizontal acceleration. The sampling frequency for the ground speed parameters includes 1Hz, and the sampling frequency for the longitudinal acceleration parameters includes 8Hz. The ground speed parameters are derived from the ground speed values ​​of the Inertial Reference System (IRS), and are determined by comparing and verifying the ground speed values ​​recorded from different data sources (GPS, IRS, and Display Management Computer (DMC)). The initial velocity vn is taken as the velocity value of the 1Hz ground velocity data point, in units of knots; Between adjacent 1Hz ground velocity data points, interpolation is performed using the above interpolation formula, Δt=0.125 seconds, or 1 / 8 seconds, to divide the 1-second interval into 8 interpolation points. By increasing the sampling frequency of the ground speed parameter from 1Hz to 8Hz, the ground speed and the vertical overload parameter at 8Hz are perfectly aligned in time. This fully utilizes the 8Hz sampling frequency of the vertical overload. Compared to the traditional 1Hz ground speed integration algorithm, which may have an error exceeding ±20 meters due to the large integration interval, increasing the sampling frequency of the ground speed parameter from 1Hz to 8Hz before integration improves the positioning accuracy by 60% (from ±20 meters to ±8 meters). Table 2 lists examples of calculating ground speed and longitudinal acceleration using the acceleration-assisted interpolation method. Table 2:

[0023] The vertical velocity is obtained by processing the vertical overload, and then the equivalent depth of the pit in the uneven section of the first runway surface is calculated according to the energy conservation model. In this embodiment, the vertical velocity is obtained by processing the vertical overload, and then the equivalent depth of the pits in the uneven section of the first runway surface is calculated according to the energy conservation model. The method is as follows: The vertical acceleration is obtained by performing gravity-based correction on the original vertical overload. , ; In the formula, The parameters are: vertical overload, 1.0 is the gravity reference value, t is the time variable, and to eliminate high-frequency noise, ... Low-pass filtering (using a Butterworth low-pass filter) is performed. The cutoff frequency fc is used to perform low-pass filtering on the vertical acceleration according to the aircraft type. The low-pass filtering process includes: for Boeing 737 series aircraft, a first cutoff frequency is used ( For the Airbus A320 family of aircraft, a second cutoff frequency is used ( For other aircraft types, the third cutoff frequency is used. After filtering, the filtered vertical acceleration at time t is obtained. ; Due to the lag in the landing gear's dynamic response, the filtered vertical acceleration is integrated with phase compensation to obtain the vertical velocity in the body coordinate system. The phase compensation integral processing includes: a first phase compensation angle for hydropneumatic strut landing gear and a second phase compensation angle for rocker arm landing gear, expressed as: ; in, The imaginary unit; The phase compensation angle is determined according to the landing gear type. For example, the phase compensation angle for a hydropneumatic shock absorber strut is 15°, and the phase compensation angle for a rocker arm landing gear is 10°. This is a phase compensation factor used to compensate for the phase lag in the landing gear's dynamic response. For at any time The filtered vertical acceleration, The dummy variable for integration (a time variable, ranging from t0 to t) represents the various time points at which the integrand takes its value; Based on aircraft pitch angle and roll angle The vertical velocity in the body coordinate system Vertical velocity converted to ground coordinates , ; Based on the peak vertical velocity in the ground coordinate system, combined with the aircraft mass and effective elastic coefficient, and according to the principle of conservation of kinetic and elastic potential energy, the equivalent depth of the pit in the uneven section of the first runway surface is calculated. To more comprehensively assess the hazards, a kinetic-potential energy conversion factor K is introduced as an intermediate variable, and a normalized roughness index is constructed. First, the kinetic-potential energy conversion factor K is calculated. The K value directly reflects the magnitude of the impact kinetic energy relative to the reference potential energy, thus representing the depth of the impact. The formula for K is: ; in, for peak value For reference depth (e.g., the statistical average of the static compression of the shock absorber strut under typical landing weight of a domestically produced narrow-body aircraft B737 / A320, preferably 0.05m), For aircraft quality, The effective elastic coefficient is matched from a preset database based on the aircraft type and runway pavement type (e.g., asphalt / concrete). The value, for example: narrow-body machines on asphalt pavement: If no specific information is available, the default value will be used. ; It also includes directly estimating the equivalent deformation depth h of the pavement based on the principle of conservation of kinetic and elastic potential energy: ; After calculating the equivalent depth of the craters within the uneven section of the first runway surface, a comprehensive evaluation is also included, specifically: Based on the equivalent depth, combined with the average vertical acceleration and the calculated ratio of the uneven length to the typical length, the comprehensive unevenness index is calculated by weighted summation. ; Where α, β, and γ are all weighting coefficients. In uneven areas Internal, the absolute value of the filtered vertical acceleration The average value, Typical length refers to the reference length of the unevenness of the corresponding type (long-wave undulation type, mixed type, short-wave impact type), and L is the calculated unevenness length; The comprehensive assessment also includes a grading system, specifically: The unevenness level is classified according to the comprehensive unevenness index, including: when the comprehensive unevenness index is less than the first level threshold (0.1), it is rated as Level 1 slight unevenness, and it is recommended to strengthen monitoring; when the comprehensive unevenness index is greater than or equal to the first level threshold and less than the second level threshold (0.3), it is rated as Level 2 moderate unevenness, and it is recommended to include it in planned inspections; when the comprehensive unevenness index is greater than or equal to the second level threshold and less than the third level threshold (0.6), it is rated as Level 3 severe unevenness, and it is recommended to prioritize repairs; when the comprehensive unevenness index is greater than or equal to the third level threshold, it is rated as Level 4 dangerous unevenness, and it is recommended to inspect it immediately and take restrictive measures. S5, Based on the length, locate the position of the pit in the uneven section of the first runway surface; In this embodiment, the method for locating the spatial position of the pit within the uneven section of the first runway surface in the runway is as follows: S5.1, Based on the high-frequency ground speed sequence, calculate the cumulative displacement of each sampling point relative to the starting point of the interval within the first runway surface unevenness interval. ( ,in, The time interval corresponding to 8Hz sampling is 0.125 seconds. (where k is the high-frequency ground velocity at the k-th sampling point, in m / s). S5.2, Based on GPS sampling data, calculate the geographic coordinates at any time within the uneven section of the first runway surface using a time interpolation method; ; ; In the formula, The latitude value at the target time t, in degrees (°), where t1 and t2 are two adjacent GPS sampling times, satisfying t1 <t<t2, The GPS latitude value at time t1. The GPS latitude value at time t2; The longitude value at the target time t is expressed in degrees (°). The GPS longitude value at time t1. The GPS longitude value at time t2; S5.3, establish a runway coordinate system with the runway entrance center point as the origin and the direction of the runway centerline as the positive direction (with the direction of the runway centerline as the positive X-axis and the direction perpendicular to the runway centerline as the Y-axis); S5.4, Combining the cumulative displacement, the true runway heading, and the geographic coordinates, calculate the coordinates of the start and end points of the first runway surface unevenness section in the runway coordinate system to obtain the spatial position of the pit in the first runway surface unevenness section within the runway. ; ; ; ; in, The coordinates of the starting point of the uneven section in the runway coordinate system are given in meters (m). The coordinates of the end point of the uneven section in the runway coordinate system are given in meters (m). This is the distance from the track entrance to the GPS starting point of this run. The true heading of the runway (unit: degrees (°), representing the direction of the runway centerline, used to project the cumulative displacement onto the X and Y axes of the runway coordinate system). The starting time in the uneven interval The corresponding cumulative displacement, in meters (m), is the distance from the GPS starting point to the beginning of the uneven section. The end time in the uneven interval The corresponding cumulative displacement, in meters (m), is the distance from the GPS starting point to the end of the uneven section. The direction cosine is used to convert the cumulative displacement from the geographic direction to the X and Y components in the runway coordinate system; This method ultimately yields the precise spatial location of the uneven area within the runway; in this example, the location result shows that the center of the uneven area is located approximately 700 meters from the entrance of runway 07. Figure 7); The steps S5.1 to S5.4 above provide a method for high-precision positioning of uneven areas detected in a single flight. As a supplementary verification or rapid macroscopic positioning method, a method based on a large amount of historical data statistics can also be used to identify and locate typical uneven areas on the runway. This method is particularly suitable for long-term monitoring and trend analysis. The location where takeoff thrust is engaged during aircraft takeoff is random, which may introduce errors in a single analysis due to changes in the baseline, as the traditional "distance from the runway threshold" positioning method is prone to errors. Using GPS latitude and longitude positioning can directly obtain absolute geographic coordinates, providing a unified baseline for multi-source data aggregation analysis. A macro-statistical analysis was conducted using runway 07 at Ganzhou Airport as an example: Runway parameters: 1 track (07 / 25 direction), 2600 meters long, 45 meters wide; Data sample: A total of 3,800 takeoff runway data points were collected (including 2,100 for Boeing aircraft and 1,700 for Airbus aircraft). Analysis method: Extract the GPS latitude and longitude coordinates of the time when the maximum value of VRTG parameter (usually corresponding to the maximum deviation of VRTG parameter) occurs during the takeoff roll phase of each flight segment, mark the data of different aircraft types separately (e.g., Boeing aircraft are marked in pink, and Airbus aircraft are marked in yellow) and overlay them; Statistical analysis revealed that, through coordinate overlay analysis, the maximum VRTG values ​​were concentrated in specific areas of the runway (e.g., Figure 7 The area shown in the red box indicates that the pavement features in this area continuously cause significant vertical impacts. Coordinate calculation: based on the electronic map of Ganzhou Airport ( Figure 8 Based on the coordinates and scale (1:10000), the geographical coordinates of the center point of this concentrated area can be estimated to be approximately: latitude 25°51′N, longitude 114°56′E. Therefore, the straight-line distance from this area to the entrance of runway 07 is estimated to be approximately 700 meters, located in the slightly rear-central part of the runway. This statistical result corroborates the results obtained from analyzing a single flight using the precise positioning method described in S5.1 to S5.4 (approximately 700 meters from the runway threshold), verifying the consistency of the method of the present invention at different analytical granularities. The statistical method can serve as an effective tool for quickly screening potential problem areas on the runway, while the high-precision positioning method provides precise spatial guidance for specific maintenance operations. By introducing a series of key technologies, including multi-resolution wavelet pattern recognition, aircraft dynamics parameter adaptation, acceleration-assisted high-frequency reconstruction, a dual-energy comprehensive assessment model, a precise positioning method based on high-frequency ground speed and GPS interpolation, and multi-segment statistical verification, this method successfully transforms massive amounts of routine flight data into a dynamic, quantitative, and localizable deep perception of runway pavement health status. This represents a fundamental shift from "periodic static detection" to "dynamic monitoring along with flight operations," significantly improving the initiative, precision, and economy of runway safety management.

[0024] Those skilled in the art will understand that, for the sake of convenience and brevity, the specific working process and related explanations of the methods described above can be found in the corresponding processes in the foregoing system embodiments, and will not be repeated here.

[0025] Those skilled in the art will recognize that the modules and method steps of the various examples described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, computer software, or a combination of both. The programs corresponding to the software modules and method steps can be placed in random access memory (RAM), main memory, read-only memory (ROM), electrically programmable ROM, electrically erasable programmable ROM, registers, hard disks, removable disks, CD-ROMs, or any other form of storage medium known in the art. To clearly illustrate the interchangeability of electronic hardware and software, the components and steps of the various examples have been generally described in terms of functionality in the foregoing description. Whether these functions are implemented in electronic hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementation should not be considered beyond the scope of the invention.

[0026] The terms “first”, “second”, etc., are used to distinguish similar objects, not to describe or indicate a specific order or sequence.

[0027] The term "comprising" or any other similar term is intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus / device that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent in such process, method, article, or apparatus / device.

[0028] The technical solution of the present invention has been described above with reference to the preferred embodiments shown in the accompanying drawings. However, it will be readily understood by those skilled in the art that the scope of protection of the present invention is obviously not limited to these specific embodiments. Without departing from the principles of the present invention, those skilled in the art can make equivalent changes or substitutions to the relevant technical features, and the technical solutions after such changes or substitutions will all fall within the scope of protection of the present invention.

Claims

1. A method for runway pavement smoothness identification based on flight data, characterized in that, The method includes the following steps: S1, obtain the QAR dataset of the target aircraft during the takeoff and landing taxiing phase of the first runway of the target airport and extract the VRTG parameters to generate a time series set X1 of VRTG parameters; the parameters extracted from the QAR dataset also include ground speed, longitudinal overload, and vertical overload; the first runway is any runway of the target airport; S2, perform first-order difference processing on X1 to obtain the difference time series X2; S3, perform anomaly detection on X2 to obtain abnormal data of VRTG parameters; Based on the waveform changes of the abnormal data of the VRTG parameters, the vibration pattern is obtained, and then the uneven area of ​​the first runway surface is obtained. S4. Calculate the horizontal acceleration using the longitudinal overload of the target aircraft, and then use acceleration-assisted interpolation to complete the frequency of the ground speed parameter of the target aircraft to the frequency of the VRTG parameter. Based on the ground speed parameters after frequency completion, the length of the pit in the uneven section of the first runway surface is obtained by ground speed integration calculation. The vertical velocity is obtained by processing the vertical overload, and then the equivalent depth of the pit in the uneven section of the first runway surface is calculated according to the energy conservation model. S5, based on the length, locate the position of the pit within the uneven section of the first runway surface.

2. The method for runway pavement smoothness identification based on flight data according to claim 1, characterized in that, The method for obtaining abnormal data of the VRTG parameter for the X2 anomaly detection is as follows: Calculate the probability density P of each difference value in X2; Each probability density P is compared with a preset probability density threshold P0. If P≤P0, the VRTG parameter of the next time point in X1 corresponding to the difference value is determined to be abnormal data, and thus an abnormal data set is obtained.

3. The method for runway pavement smoothness identification based on flight data according to claim 2, characterized in that, Based on the waveform changes of the abnormal data of the VRTG parameters, the vibration pattern is obtained, and then the unevenness range of the first runway surface is determined. The unevenness range includes three types: long-wave undulation unevenness, mixed unevenness, and short-wave impact unevenness. The method is as follows: S3.1 Calculate the degree of deviation for each abnormal data in the abnormal data set; the degree of deviation is the standardized deviation of the probability density value corresponding to the abnormal data relative to the statistical distribution of the probability density of normal samples; S3.2, Perform multi-resolution wavelet decomposition on the original VRTG parameter time series corresponding to the abnormal data set to extract low-frequency and high-frequency feature components of a set Hertz. S3.3 Calculate the proportion of the energy of the low-frequency characteristic component in the total characteristic energy. When the proportion is greater than the first threshold, it is determined to be a long-wave undulation type unevenness. When the proportion is between the first threshold and the second threshold, it is determined to be a mixed type unevenness. When the proportion is less than the second threshold, it is determined to be a short-wave impact type unevenness. And determine the initial time length according to the determined unevenness type. S3.4, taking the abnormal data with the largest deviation as the center, extend to both sides by half of the initial time length to form the initial abnormal interval; S3.5 Calculate the comprehensive anomaly index for each time point within the initial anomaly interval. The comprehensive anomaly index is a weighted combination value based on the degree of deviation, the rate of change of VRTG parameters, and local volatility. The weight coefficients of the weighted combination are adaptively adjusted according to the unevenness type. S3.6, start from the boundary of the initial abnormal interval and iteratively expand outward. Determine whether to include adjacent time points in the uneven interval based on the comparison between the comprehensive abnormal index and the preset threshold. Stop expanding when the comprehensive abnormal index of multiple consecutive time points is lower than the preset threshold. The preset threshold is dynamically determined according to the statistical characteristics of the comprehensive abnormal index in the initial abnormal interval. S3.7 determines the final uneven interval after iterative expansion of the time boundary and outputs the corresponding unevenness type.

4. The runway pavement smoothness identification method based on flight data according to claim 1, characterized in that, The frequency of the ground speed parameter of the target aircraft is supplemented to the frequency of the VRTG parameter using an acceleration-assisted interpolation method. The method is as follows: Based on the initial velocity and the current horizontal acceleration of the ground speed parameter, interpolation calculations are performed according to a preset time interval. The sampling interval of the target aircraft is divided into multiple interpolation points for interpolation completion, so that the sampling frequency of the ground speed parameter reaches the frequency of the VRTG parameter.

5. The method for runway pavement smoothness identification based on flight data according to claim 1, characterized in that, The vertical velocity is obtained by processing the vertical overload, and then the equivalent depth of the pits in the uneven section of the first runway surface is calculated according to the energy conservation model. The method is as follows: The vertical acceleration is obtained by performing gravity-based correction on the original vertical overload. The vertical acceleration is subjected to low-pass filtering according to the aircraft type. The low-pass filtering includes: a first cutoff frequency for Boeing 737 series aircraft, a second cutoff frequency for Airbus A320 series aircraft, and a third cutoff frequency for other aircraft types. The filtered vertical acceleration is subjected to phase compensation integration processing to obtain the vertical velocity in the body coordinate system. The phase compensation integration processing includes: using a first phase compensation angle for oil-gas type shock absorber strut landing gear and a second phase compensation angle for rocker arm type landing gear. Based on the aircraft's pitch and roll angles, the vertical velocity in the body coordinate system is converted into the vertical velocity in the ground coordinate system; Based on the peak vertical velocity in the ground coordinate system, combined with the aircraft mass and effective elastic coefficient, and according to the principle of conservation of kinetic and elastic potential energy, the equivalent depth of the pits in the uneven section of the first runway surface is calculated.

6. The runway pavement smoothness identification method based on flight data according to claim 5, characterized in that, After calculating the equivalent depth of the craters within the uneven section of the first runway surface, a comprehensive evaluation is also included, specifically: Based on the equivalent depth, combined with the average vertical acceleration and the calculated ratio of the uneven length to the typical length, the comprehensive unevenness index is calculated by weighted summation.

7. The runway pavement smoothness identification method based on flight data according to claim 6, characterized in that, The comprehensive assessment also includes a grading system, specifically: The unevenness level is classified according to the comprehensive unevenness index, including: when the comprehensive unevenness index is less than the first level threshold, it is rated as Level 1 slight unevenness; when the comprehensive unevenness index is greater than or equal to the first level threshold and less than the second level threshold, it is rated as Level 2 moderate unevenness; when the comprehensive unevenness index is greater than or equal to the second level threshold and less than the third level threshold, it is rated as Level 3 severe unevenness; and when the comprehensive unevenness index is greater than or equal to the third level threshold, it is rated as Level 4 dangerous unevenness.

8. The method for runway pavement smoothness identification based on flight data according to claim 1, characterized in that, The ground velocity parameters are completed to a sampling frequency consistent with the VRTG parameters using an acceleration-assisted interpolation method. The method is as follows: Based on the initial velocity and the current horizontal acceleration of the ground speed parameters, interpolation calculations are performed according to a preset time interval. The sampling interval of the target aircraft is divided into multiple interpolation points for interpolation completion, so that the sampling frequency of the ground speed parameters reaches the frequency of the VRTG parameters, thus obtaining a high-frequency ground speed sequence.

9. The method for runway pavement smoothness identification based on flight data according to claim 1, characterized in that, The method for performing first-order differencing on the time series set X1 to obtain the differencing time series X2 is as follows: Calculate the first difference of each element in the time series set X1 to obtain the differenced time series X2.

10. The method for runway pavement smoothness identification based on flight data according to claim 8, characterized in that, The method for locating the spatial position of the pit within the uneven section of the first runway surface is as follows: S5.1, Based on the high-frequency ground speed sequence, calculate the cumulative displacement of each sampling point relative to the starting point of the interval within the uneven section of the first runway surface; S5.2, Based on GPS sampling data, calculate the geographic coordinates at any time within the uneven section of the first runway surface using a time interpolation method; S5.3, establish a runway coordinate system with the runway entrance center point as the origin and the runway centerline direction as the positive direction; S5.4, Combining the cumulative displacement, the true runway heading, and the geographic coordinates, calculate the coordinates of the start and end points of the first runway surface unevenness section in the runway coordinate system to obtain the spatial location of the pits in the first runway surface unevenness section within the runway.

Images