Spatial Continuous Visibility Prediction Method and Device Based on Multimodal Data Analysis

Abstract

The present invention provides a method and device for spatially continuous visibility prediction based on multi-modal data analysis. The method includes: obtaining multi-modal observation data of a target area, performing image feature encoding on the image observation data therein to generate atmospheric optical thickness parameters, and performing path feature encoding on the lidar scan data therein to generate extinction coefficient distribution parameters; inputting the atmospheric optical thickness parameters and the extinction coefficient distribution parameters into a multi-modal fusion model to perform horizontal optical thickness correction on the atmospheric optical thickness parameters to generate horizontal optical thickness correction parameters, and performing spatial continuity fusion on the horizontal optical thickness correction parameters and the extinction coefficient distribution parameters to generate the spatial continuous visibility distribution characteristics of the target area, and then converting the spatial continuous visibility distribution characteristics into a visibility prediction result in the polar coordinate system, and outputting a coordinate distribution map of the visibility prediction result. By means of the present invention, the accuracy of visibility prediction can be improved.

Description

Technical Field

[0001] The present invention relates to the field of data processing, particularly to the field of meteorological data processing, and specifically to a method and device for predicting spatially continuous visibility based on multi-modal data analysis. Background Art

[0002] Visibility prediction is an important technical means in the fields of meteorological monitoring and environmental perception, and it plays a crucial role in scenarios such as aviation safety, transportation, and disaster warning. Traditional visibility prediction methods usually rely on the observation and analysis of a single data source. For example, the atmospheric optical parameters are retrieved from visible light images, or the aerosol extinction coefficient is estimated based on the horizontal scan data of lidar. However, such methods have significant limitations in practical applications and are difficult to achieve spatially continuous and high-precision visibility prediction. Currently, lidar observes visibility by horizontally scanning at a fixed angle to obtain the path visibility distribution. In order to obtain spatially continuous visibility, generally, a set of horizontally scanned data is directly interpolated to obtain the spatial visibility. However, this method has a large error in direct interpolation at locations far from the lidar visibility observation instrument, where the distance between two observations is relatively large. Summary of the Invention

[0003] The present invention provides a method and device for predicting spatially continuous visibility based on multi-modal data analysis.

[0004] According to one aspect of the present invention, there is provided a method for predicting spatially continuous visibility based on multi-modal data analysis. The method includes: acquiring multi-modal observation data of a target area, where the multi-modal observation data includes image observation data and lidar scan data, and the image observation data includes visible light images of the target area and corresponding shooting time and spatial position information, and the lidar scan data includes aerosol extinction coefficient distribution data on the horizontal scan path of the target area; performing image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and performing path feature encoding on the lidar scan data to generate extinction coefficient distribution parameters; inputting the atmospheric optical thickness parameters and the extinction coefficient distribution parameters into a multi-modal fusion model, and performing horizontal optical thickness correction on the atmospheric optical thickness parameters through the multi-modal fusion model to generate horizontal optical thickness correction parameters; performing spatial continuity fusion based on the horizontal optical thickness correction parameters and the extinction coefficient distribution parameters to generate spatially continuous visibility distribution features of the target area; converting the spatially continuous visibility distribution features into a visibility prediction result in the polar coordinate system, and outputting a coordinate distribution map of the visibility prediction result.

[0005] According to another aspect of the present invention, there is provided a spatial continuous visibility prediction device based on multi-modal data analysis, comprising: a data acquisition module for acquiring multi-modal observation data of a target area, the multi-modal observation data including image observation data and lidar scan data, wherein the image observation data includes visible light images of the target area and corresponding shooting time and spatial position information, and the lidar scan data includes aerosol extinction coefficient distribution data on the horizontal scan path of the target area; a feature encoding module for performing image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and performing path feature encoding on the lidar scan data to generate extinction coefficient distribution parameters; a parameter correction module for inputting the atmospheric optical thickness parameters and the extinction coefficient distribution parameters into a multi-modal fusion model, and performing horizontal optical thickness correction on the atmospheric optical thickness parameters through the multi-modal fusion model to generate horizontal optical thickness correction parameters; a parameter fusion module for performing spatial continuity fusion based on the horizontal optical thickness correction parameters and the extinction coefficient distribution parameters to generate spatial continuous visibility distribution characteristics of the target area; and a visibility prediction module for converting the spatial continuous visibility distribution characteristics into a visibility prediction result in the polar coordinate system and outputting a coordinate distribution map of the visibility prediction result.

[0006] The present invention has at least the following beneficial effects: The spatial continuous visibility prediction method based on multi-modal data analysis provided by the present invention fuses image observation data and lidar scan data, extracts the atmospheric optical thickness parameter and the aerosol extinction coefficient distribution parameter respectively, and performs horizontal optical thickness correction and spatial continuity fusion on the two types of parameters based on the multi-modal fusion model, and finally generates the visibility prediction result in the polar coordinate system. Based on this, the image observation data can reflect the atmospheric scattering characteristics through the optical attenuation gradient of the visible light image, and the lidar scan data can provide the radial distance details of the aerosol concentration through the extinction coefficient distribution on the horizontal scan path. The collaborative processing of the two can make full use of the advantages of the image in the horizontal angular resolution and the lidar in the radial distance observation, and make up for the problem of insufficient spatial coverage of a single data source. Moreover, when performing horizontal correction on the atmospheric optical thickness parameter based on the multi-modal fusion model, the weight distribution of the optical thickness is dynamically optimized through the spatial correlation of the lidar extinction coefficient distribution, which can suppress the errors caused by light changes or local occlusion in the image data and improve the accuracy of physical parameter coupling. In addition, the corrected optical thickness parameter and the extinction coefficient distribution parameter are fused in space continuity in the polar coordinate system to generate a seamless visibility distribution feature covering both the angle and distance dimensions, avoiding the information distortion problem in the traditional Cartesian coordinate system conversion, and making the prediction result more conform to the original observation characteristics of the lidar and the camera. Through the output of the visibility distribution map in the polar coordinate system, the spatial continuous visibility change in the target area within the horizontal scan path and the image observation angle range can be directly reflected, significantly improving the interpretability of the prediction result and the decision-making support efficiency in actual meteorological services. Description of the Drawings

[0007] Figure 1 Fig. shows a schematic diagram of an application scenario of a spatial continuous visibility prediction method based on multi-modal data analysis according to an embodiment of the present invention.

[0008] Figure 2 Fig. shows a flowchart of a spatial continuous visibility prediction method based on multi-modal data analysis according to an embodiment of the present invention.

[0009] Figure 3 Fig. shows a schematic diagram of the functional module architecture of a spatial continuous visibility prediction device based on multi-modal data analysis according to an embodiment of the present invention. Detailed Embodiments

[0010] Figure 1A schematic diagram of an application scenario provided according to an embodiment of the present invention is shown. The application scenario includes one or more data acquisition devices 101, a data processing device 120, and one or more networks 110 that couple the one or more data acquisition devices 101 to the data processing device 120. The data acquisition device 101 at least includes an image acquisition device (such as a pan-tilt camera) and a lidar.

[0011] In Figure 1 the configuration shown, the data processing device 120 may include one or more components that implement the functions performed by the data processing device 120. These components may include software components, hardware components, or a combination thereof that may be executed by one or more processors. It should be understood that various different system configurations are possible, which may be different from the application scenario. Therefore, Figure 1 is an example of a system for implementing the various methods described herein and is not intended to be limiting. The above application scenario may also include one or more databases 130. In some embodiments, these databases may be used to store data and other information. For example, one or more of the databases 130 may be used to store meteorological observation data, such as image data and lidar data. The database 130 may reside at various locations. For example, the database used by the data processing device 120 may be local to the data processing device 120, or may be remote from the data processing device 120 and may communicate with the data processing device 120 via a network-based or dedicated connection.

[0012] Please refer to Figure 2 , the spatial continuous visibility prediction method based on multimodal data analysis provided by the embodiment of the present invention may specifically include the following steps:

[0013] Step S100: Obtain multimodal observation data of the target area. The multimodal observation data includes image observation data and lidar scan data. The image observation data includes visible light images of the target area and corresponding shooting time and spatial position information, and the lidar scan data includes aerosol extinction coefficient distribution data on the horizontal scan path of the target area.

[0014] The target area is, for example, the spatial range where visibility prediction is required. This area covers continuous observation points within a specific geographical range, such as the continuous space where the urban road network, the periphery of the airport runway, or meteorological monitoring stations are located. The multi-modal observation data is a set of heterogeneous data collected through different sensors or observation means, specifically including two types in this step: image observation data and lidar scan data. The image observation data is captured by a visible light imaging device (such as a pan-tilt camera), containing two-dimensional image information of the target area in the visible light band, along with the shooting time and spatial position information. The shooting time is used to identify the precise moment of image acquisition (e.g., year, month, day, hour, minute, and second in Coordinated Universal Time format), and the spatial position information is used to record the longitude and latitude coordinates of the imaging device and the shooting angle parameters, ensuring that the image data can be mapped to the geographical space coordinate system. The visible light image is a multi-channel pixel matrix, and each pixel point records the light intensity information of the corresponding position in the target area, which can be used to analyze the atmospheric transmission characteristics. The lidar scan data is generated by a lidar device through actively emitting laser pulses and receiving the backscattered signals, containing the aerosol extinction coefficient distribution data along the horizontal scan path in the target area. The horizontal scan path refers to the spatial trajectory formed by the lidar rotating periodically in the horizontal direction at a fixed altitude. The aerosol extinction coefficient distribution data characterizes the degree of energy attenuation caused by the absorption and scattering of aerosol particles in the propagation path, and its value has a negative correlation with the atmospheric visibility. Specifically, the spatial position information of the image observation data needs to be spatially registered with the horizontal scan path of the lidar scan data to ensure that both cover the same geographical range and observation direction. For example, when the lidar performs a horizontal scan along the due north direction, the shooting angle of the visible light image needs to be adjusted to the due north direction, and the image pixels are aligned with the azimuth angle of the lidar scan path based on the longitude and latitude coordinates in the spatial position information. Further, the aerosol extinction coefficient distribution data in the lidar scan data is discretely sampled according to the distance resolution. For example, an extinction coefficient value is recorded every 10 meters, forming a continuous extinction coefficient sequence along the horizontal scan path, while the visible light image in the image observation data covers the two-dimensional space of the same area according to a preset resolution (such as 1920×1080 pixels). By synchronously acquiring the two types of data, a multi-modal observation data set of the target area can be constructed, providing spatially aligned inputs for subsequent feature extraction and fusion.

[0015] Step S200: Perform image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and perform path feature encoding on the lidar scan data to generate extinction coefficient distribution parameters.

[0016] Image feature encoding is the process of extracting feature parameters related to atmospheric optical properties from visible light images through computer vision algorithms. The atmospheric optical thickness parameter is used to quantify the attenuation degree of the atmosphere to visible light, and its value is determined by the concentration and distribution of components such as aerosols and water vapor. Specifically, it is defined as the total extinction integral value when light passes through the atmosphere along a vertical path. In this step, the generation of the atmospheric optical thickness parameter needs to combine the brightness distribution and optical attenuation gradient in the sky region of the visible light image. For example, by segmenting the sky region in the visible light image (such as identifying the blue pixels in the top region of the image based on color thresholds or deep learning models), the mean pixel brightness in the sky region at different azimuth angles (such as taking due north as 0 degrees and dividing a fan-shaped region every 15 degrees clockwise) is extracted to form a pixel brightness distribution sequence. Based on the change trend of the brightness means between adjacent azimuth angles in this sequence (such as the brightness decreasing with the increase of the azimuth angle), the optical attenuation gradient distribution in the horizontal direction can be calculated, reflecting the spatial difference in light intensity caused by atmospheric inhomogeneity. Further, the optical attenuation gradient distribution is input into a pre-trained atmospheric scattering correction network. This network models the non-linear relationship between the optical attenuation gradient and atmospheric scattering parameters (such as Rayleigh scattering coefficient, Mie scattering coefficient) through a multi-layer perceptron structure and outputs the initial atmospheric optical thickness distribution. To eliminate the influence of the solar elevation angle on the optical thickness calculation, the solar elevation angle parameter corresponding to the shooting time needs to be introduced (such as calculating the elevation angle of the sun relative to the horizon at the shooting moment through astronomical algorithms), and the time correction module in the atmospheric scattering correction network is used to perform angle compensation on the initial atmospheric optical thickness distribution to generate the corrected optical thickness distribution, aligning it with the actual solar azimuth. Finally, the corrected optical thickness distribution is spatially projected according to the fan-shaped observation sub-regions of the target area (such as the fan-shaped regions divided every 15 degrees), and local atmospheric optical thickness parameters are assigned to each fan-shaped sub-region. At the same time, the path feature encoding of lidar scan data refers to extracting statistical features related to visibility from the aerosol extinction coefficient distribution data of the horizontal scan path. Specifically, the horizontal scan path is segmented into multiple radial distance segments according to a preset distance resolution (such as every 50 meters). For each distance segment, the mean extinction coefficient of the continuous extinction coefficient sampling points inside it is calculated, and the extinction integral contribution value is generated based on the product of the mean extinction coefficient and the length of the distance segment, reflecting the cumulative impact of this distance segment on the overall visibility. The extinction integral contribution values of all distance segments are arranged in spatial order to form an extinction integral value sequence, and then the extinction integral values are mapped to extinction coefficient distribution parameters matching the fan-shaped observation sub-regions through a visibility conversion relationship (such as the Koschmieder formula). Through the above processing, the image observation data and lidar scan data are respectively encoded into spatially aligned atmospheric optical thickness parameters and extinction coefficient distribution parameters, providing standardized inputs for multi-modal fusion.

[0017] Step S300: Input the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into the multi-modal fusion model. Through the multi-modal fusion model, perform horizontal optical thickness correction on the atmospheric optical thickness parameter to generate a horizontal optical thickness correction parameter.

[0018] The multi-modal fusion model is, for example, a feature alignment and cross-modal association model designed based on a deep learning architecture, used to integrate feature parameters from different modalities and eliminate the inconsistencies between modalities. The horizontal optical thickness correction refers to correcting the deviation of the atmospheric optical thickness parameter derived from image observation data in the horizontal angle according to the extinction coefficient distribution parameter provided by lidar scan data. Specifically, although the atmospheric optical thickness parameter can reflect the extinction characteristics of the atmosphere in the vertical direction, its horizontal resolution is limited by the shooting angle of the visible light image and the accuracy of the optical attenuation gradient calculation, while the extinction coefficient distribution parameter provided by lidar scan data has a high-precision horizontal spatial resolution. Therefore, the feature alignment layer of the multi-modal fusion model first converts the atmospheric optical thickness parameter into a first feature vector in the horizontal angle dimension (for example, expanding the parameters of each 15-degree sector sub-region into a one-dimensional vector), and at the same time converts the extinction coefficient distribution parameter into a second feature vector in the distance dimension (for example, expanding the parameters of each 50-meter distance segment into another one-dimensional vector). Calculate the association weight between the first feature vector and the second feature vector through the cross-attention mechanism. For example, for the atmospheric optical thickness parameter of a certain sector sub-region, the model automatically focuses on the extinction coefficient distribution parameter within the lidar distance segment corresponding to the sector angle, and performs weighted fusion on the first feature vector based on the association weight to generate a fused feature vector. After the fused feature vector is mapped through the fully connected layer, a horizontal optical thickness correction parameter is output. This parameter, while retaining the vertical extinction information of the atmosphere, integrates the high-resolution extinction coefficient data of the lidar horizontal scan path, thus significantly improving the accuracy of the optical thickness estimation in the horizontal direction. For example, if the initial estimated value of the atmospheric optical thickness parameter of a certain sector sub-region is too high due to cloud cover, and the lidar detects a low extinction coefficient value (indicating high visibility) within the distance segment corresponding to this angle, the model reduces the contribution of the optical thickness parameter of this sector sub-region through the attention weight to generate a more realistic horizontal optical thickness correction parameter.

[0019] Step S400: Perform spatial continuity fusion based on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter to generate the spatial continuous visibility distribution feature of the target area.

[0020] Spatial continuity fusion is to convert discrete feature parameters into a continuous visibility distribution covering the entire space of the target area through interpolation and smoothing. Although the horizontal optical thickness correction parameter has been calibrated by a multi-modal fusion model, its spatial resolution is still limited by the division of fan-shaped sub-regions (such as one fan every 15 degrees), while the spatial resolution of the extinction coefficient distribution parameter is limited by the division of lidar distance segments (such as one distance segment every 50 meters). To achieve spatially continuous visibility estimation, the two types of parameters need to be linearly interpolated on an angle-distance grid. Specifically, the horizontal optical thickness correction parameter is divided into finer-grained angular sub-parameters according to a preset angular resolution (such as every 1 degree), and at the same time, the extinction coefficient distribution parameter is divided into distance sub-parameters according to a preset distance resolution (such as every 10 meters); bilinear interpolation is performed on each angular sub-parameter and the corresponding distance sub-parameter to generate the visibility interpolation result within the continuous angle-distance grid. For example, at a position with an azimuth angle of 30.5 degrees and a distance of 55 meters, its visibility value is calculated by the weighted average of the angular sub-parameters of 30 degrees and 31 degrees and the distance sub-parameters of 50 meters and 60 meters. Subsequently, the visibility interpolation result is input into the spatial continuity optimization network, which smooths the visibility values of adjacent grids through convolutional layers and pooling layers to eliminate local mutation noise caused by interpolation. For example, if the visibility value of a grid cell is significantly higher than that of its surrounding cells, the network reduces the impact of this outlier through gradient optimization to generate a spatially continuous visibility distribution feature. This feature covers the entire space of the target area in the form of a high-resolution grid and can accurately characterize the continuous trend of visibility changing with azimuth angle and distance.

[0021] Step S500: Convert the spatially continuous visibility distribution feature into a visibility prediction result in the polar coordinate system and output the coordinate distribution map of the visibility prediction result.

[0022] The visibility prediction result in the polar coordinate system refers to a two-dimensional matrix with angle-distance as the coordinate axes, where each matrix cell stores the visibility value corresponding to the azimuth angle and distance. The coordinate distribution map visualizes this matrix as a polar coordinate heat map to facilitate the intuitive display of the visibility distribution in all directions of the target area. Specifically, based on the maximum radial distance of the spatially continuous visibility distribution feature (such as the maximum coverage distance of 5 kilometers scanned by the lidar) and the preset angular interval of the fan-shaped observation sub-region (such as 15 degrees), polar coordinate grid division parameters are generated, including the angle axis division sequence (such as 0 degrees, 15 degrees, 30 degrees... 345 degrees) and the distance axis division sequence (such as 0 meters, 100 meters, 200 meters... 5000 meters). Convert the original Cartesian coordinates (such as x = 3000 meters, y = 4000 meters) of each visibility value in the spatially continuous visibility distribution feature into the angular component (such as arctan(y / x) = 53.13 degrees) and distance component (such as √(x 2 +y 2) = 5000 m), and map it to the corresponding grid cells by dividing the sequence along the angular axis and the distance axis. For the blank cells not covered in the initial polar coordinate visibility distribution mapping table (such as a grid cell without an original visibility value falling into it), the two-way adjacent interpolation method is used to extract the visibility values in the adjacent angular directions (such as the previous 15 degrees and the next 15 degrees) and the distance directions (such as the previous 100 m and the next 100 m) of this cell, and a weighted average is taken to generate the interpolated visibility value. Finally, the visibility values of all grid cells are rearranged according to the sequences of the angular axis and the distance axis to form a visibility prediction result matrix covering the complete polar coordinate range, and it is rendered as a coordinate distribution map through a visualization tool. Different colors or contour lines in the map represent the high and low visibility values. For example, the red area represents low visibility (less than 1 km), and the green area represents high visibility (greater than 10 km). This coordinate distribution map can be directly output to the meteorological warning system or the traffic management platform to provide data support for real-time decision-making.

[0023] As an implementation manner, in step S100, obtaining multi-modal observation data of the target area may specifically include:

[0024] Step S110: Extract the shooting timestamp of the image observation data and the scanning timestamp of the lidar scanning data, and mark the image observation data and the lidar scanning data with a timestamp difference less than a preset synchronization threshold as an associated data set with time alignment.

[0025] The shooting timestamp of the image observation data is the precise time identifier when the visible light image acquisition device records data (such as "2023-08-20 10:00:00.000" in Coordinated Universal Time format), and the scanning timestamp of the lidar scanning data is the starting moment when the lidar completes a single horizontal scan (such as "2023-08-20 10:00:02.500"). The preset synchronization threshold is used to define the time alignment tolerance range of the two types of data. For example, it is set to 3 seconds. If the absolute value of the difference between the shooting timestamp of the image observation data and the scanning timestamp of the lidar scanning data is less than 3 seconds, it is considered that the two are synchronized in the time dimension and cross-modal data association can be performed. Specifically, traverse the timestamps of all image observation data and lidar scanning data, calculate the timestamp difference of each pair of data, and mark the data pairs that meet the threshold conditions as an associated data set with time alignment. For example, if the shooting time of a visible light image is 10:00:01 and the starting time of the lidar scan is 10:00:03, and the timestamp difference is 2 seconds (less than 3 seconds), then this image and the scan data are classified into the same associated data set; conversely, if the timestamp difference is 4 seconds, it is excluded. Through this step, multi-modal data with time synchronization is screened out to ensure the physical consistency of subsequent feature fusion and avoid spatial feature mismatch caused by time misalignment.

[0026] Step S120: Based on the horizontal scanning angle range of the lidar scan data in the associated data set, divide the target area into multiple fan-shaped observation sub-areas at a preset angle interval, and assign a corresponding radial distance interval to each fan-shaped observation sub-area.

[0027] The horizontal scanning angle range of the lidar scan data is the azimuth angle range covered by the lidar during a single scan in the horizontal plane. For example, a complete circular scan from 0 degrees (due north) to 360 degrees (equivalent to 0 degrees). The preset angle interval is used to divide the horizontal scanning angle range into equally wide fan-shaped areas. For example, each fan-shaped observation sub-area is divided every 15 degrees, forming 24 fan-shaped sub-areas (360 / 15 = 24). Each fan-shaped observation sub-area is defined by its starting angle and ending angle. For example, the first sub-area is from 0 degrees to 15 degrees, the second is from 15 degrees to 30 degrees, and so on. The radial distance interval refers to the distance coverage range of each fan-shaped observation sub-area in the radial direction. Usually, it is divided according to the maximum effective detection distance of the lidar (such as 5000 meters) and the preset distance segmentation rule (such as every 1000 meters as an interval). For example, for the fan-shaped observation sub-area from 0 degrees to 15 degrees, its radial distance interval can be divided into multiple consecutive intervals such as 0 - 1000 meters, 1000 - 2000 meters, 2000 - 3000 meters, etc. Each interval corresponds to a section of the radial distance on the lidar scan path. Through this step, the target area is discretized into multiple angle-distance grid cells, each cell being uniquely determined by the fan-shaped angle range and the radial distance interval, providing a structured framework for subsequent multi-modal data spatial binding.

[0028] Step S130: According to the spatial position information of the image observation data in the associated data set, extract the pixel area in the visible light image that overlaps with the angle range of the fan-shaped observation sub-area as a local observation image block, and spatially bind the local observation image block to the corresponding radial distance interval.

[0029] The spatial location information of the image observation data includes the longitude and latitude coordinates of the imaging device, the shooting azimuth angle (such as 0 degrees when the lens faces due north), and the pitch angle parameter. Based on this information, each pixel in the visible light image can be mapped to a spatial position in the geographic coordinate system. For example, if the imaging device is located at 116.4 degrees longitude and 39.9 degrees latitude, the shooting azimuth angle is 90 degrees (due east), and the pitch angle is the horizontal viewing angle (0 degrees), then the center pixel of the image corresponds to the due east direction, the left pixels correspond to the northward direction, and the right pixels correspond to the southward direction. According to the angular range of the fan-shaped observation sub-region (such as 30 degrees to 45 degrees), calculate the pixel region in the visible light image that overlaps with this angular range: if the shooting azimuth angle of the imaging device is 90 degrees, the actual geographic azimuth angle corresponding to the fan-shaped observation sub-region from 30 degrees to 45 degrees is the imaging device azimuth angle (90 degrees) plus the fan-shaped angle (30 - 45 degrees), that is, actually covering 120 degrees to 135 degrees (90 + 30 to 90 + 45). Map the pixel azimuth angle of the visible light image to the geographic azimuth angle through coordinate transformation, extract the pixel region that overlaps with this fan-shaped sub-region (such as a specific rectangular region in the upper right of the image), and generate a local observation image block. Further, spatially bind the local observation image block to the corresponding radial distance interval (such as 2000 - 3000 meters), that is, label the spatial unit corresponding to this image block as "30 - 45 degree fan, 2000 - 3000 meter distance interval" to ensure strict alignment of the image data and the lidar scan data in the spatial dimension. For example, if a local observation image block is bound to the "45 - 60 degrees, 1000 - 2000 meters" unit, then all its pixels represent the visible light information within this fan angle and distance interval.

[0030] Step S140: Perform distance interval filtering on the lidar scan data within each fan-shaped observation sub-region, retain the aerosol extinction coefficient distribution data within the radial distance interval, and generate an extinction coefficient subset corresponding one-to-one with the local observation image block.

[0031] Range interval filtering is to select the extinction coefficient data belonging to a specific radial range from lidar scan data. For example, for the "0-15 degree sector, 1000-2000 meters" unit, all aerosol extinction coefficient sampling points within the azimuth angle range of 0-15 degrees and the distance range of 1000-2000 meters on the lidar scan path need to be extracted. Specifically, the lidar scan data records the sequence of extinction coefficients along the radial path at a preset distance resolution (such as one sampling point per 10 meters), and each sampling point contains the azimuth angle, distance value, and extinction coefficient value. For a certain sector sub-region (such as 15-30 degrees) and its bound radial distance interval (such as 500-1500 meters), traverse all the scan data within the sector angle range, and filter out the extinction coefficient sampling points with distance values between 500 meters and 1500 meters to form a subset of extinction coefficients. The data in this subset is sorted in ascending order of distance, for example, it contains the extinction coefficient values at 500 meters, 510 meters, 520 meters... 1500 meters, and forms a one-to-one association with the local observation image block (such as the "15-30 degrees, 500-1500 meters" image block) bound to the corresponding sector sub-region. Through this step, the lidar data is segmented into subsets of extinction coefficients that are spatially aligned with the image blocks, providing input for subsequent cross-modal feature fusion.

[0032] Step S150: Arrange the local observation image blocks and the subsets of extinction coefficients in the angular order of the sector observation sub-regions to generate a spatially aligned multi-modal observation data sequence.

[0033] The angular order is the order arranged from the starting angle of the sector observation sub-region from small to large. For example, 0-15 degrees is the first sub-region, 15-30 degrees is the second sub-region, until 345-360 degrees is the last sub-region. For each sector sub-region, extract its bound local observation image block and the subset of extinction coefficients, and arrange them in sequence as a data sequence in the angular order. For example, if the target area is divided into 24 sector sub-regions, the data sequence contains 24 elements, and each element contains a pair of "local observation image block + subset of extinction coefficients". Further, the data within each sub-region is further subdivided according to the radial distance interval. For example, for the "0-15 degrees" sector sub-region, the image blocks and the subsets of extinction coefficients corresponding to distance intervals such as 0-1000 meters and 1000-2000 meters are arranged in sequence. Through this step, the multi-modal observation data sequence maintains strict spatial alignment in both the angular and distance dimensions, forming a structured input for subsequent model processing. For example, the first element of the data sequence is the "0-15 degrees, 0-1000 meters" image block and its corresponding subset of extinction coefficients for 0-1000 meters, the second element is the "0-15 degrees, 1000-2000 meters" image block and the corresponding subset of extinction coefficients, and the subsequent elements are arranged in sequence according to the angle and distance. This sequence can be directly input into the multi-modal fusion model to ensure that the model can synchronously analyze the image and lidar features within the same spatial unit.

[0034] As an implementation, in step S200, image feature encoding is performed on the image observation data to generate the atmospheric optical thickness parameter, which may specifically include:

[0035] Step S210: Extract the pixel brightness distribution sequence in the sky region of the visible light image, where the pixel brightness distribution sequence includes the brightness means of the sky region at different azimuth angles.

[0036] The sky region in the visible light image is, for example, the part of the image that is not blocked by ground objects and mainly contains atmospheric scattered light, usually located at the top of the image or a set of pixels within a specific azimuth angle range. The pixel brightness distribution sequence is generated by dividing the sky region by azimuth angle and calculating the average brightness value of each partition, where the azimuth angle division is based on a preset angular interval (such as every 15 degrees) of the fan-shaped observation sub-region of the target area. Specifically, for the angular range (such as 30 degrees to 45 degrees) corresponding to each fan-shaped observation sub-region, the sky region pixels overlapping with this angular range are extracted from the visible light image, and the brightness mean of all pixels in this region (such as the average value after converting the RGB three channels to grayscale values) is calculated to form a brightness mean sequence indexed by azimuth angle. For example, if the sky region of the target area is divided into 24 fan-shaped partitions at 15-degree intervals from 0 degrees to 360 degrees, the pixel brightness distribution sequence contains 24 brightness mean data points, and each data point corresponds to the average brightness of a fan-shaped partition. The level of the brightness mean reflects the atmospheric scattering intensity. A low brightness mean indicates optical attenuation caused by a high aerosol concentration, while a high brightness mean corresponds to a cleaner atmosphere. Through this step, the complex brightness information in the visible light image is encoded into a quantization sequence in the azimuth angle dimension, providing input for subsequent horizontal direction optical attenuation analysis.

[0037] Step S220: Generate the optical attenuation gradient distribution of the visible light image in the horizontal direction according to the change trend of the brightness means of adjacent azimuth angles in the pixel brightness distribution sequence.

[0038] The optical attenuation gradient distribution is used to characterize the spatial variation rate of the atmospheric extinction characteristics in the horizontal direction, and its value is calculated from the difference in the average brightness of adjacent azimuth angles. Specifically, for consecutive azimuth angle partitions in the pixel brightness distribution sequence (such as 0-15 degrees, 15-30 degrees, etc.), calculate the difference in the average brightness of adjacent partitions, and divide this difference by the azimuth angle interval (such as 15 degrees) to obtain the optical attenuation gradient value in the horizontal direction. For example, if the average brightness of the 0-15 degree partition is 200 (gray value range 0-255) and the average brightness of the 15-30 degree partition is 180, the difference between the two is -20, and the corresponding optical attenuation gradient is -20 / 15 ≈ -1.33 (unit: gray value / degree). A negative gradient value indicates that the brightness decreases with the increase of the azimuth angle, that is, the atmospheric extinction effect increases; a positive gradient value indicates that the brightness increases with the increase of the azimuth angle, that is, the atmospheric transparency increases. By traversing all adjacent azimuth angle partitions, an optical attenuation gradient distribution sequence covering the entire horizontal scan range is generated, for example, containing 23 gradient values (24 azimuth angle partitions correspond to 23 adjacent gradients). This distribution can reveal the local visibility differences caused by atmospheric inhomogeneity. For example, the appearance of continuous negative gradients within a certain azimuth angle range indicates the presence of an aerosol aggregation area in that direction.

[0039] Step S230: Input the optical attenuation gradient distribution into a pre-trained atmospheric scattering correction network, and map the non-linear relationship between the optical attenuation gradient distribution and the atmospheric scattering parameters through the multi-layer perception structure in the atmospheric scattering correction network to generate an initial atmospheric optical thickness distribution.

[0040] The atmospheric scattering correction network is a regression model based on deep learning. Its input layer receives the optical attenuation gradient distribution sequence, and the hidden layer models the non-linear relationship between the gradient distribution and the atmospheric scattering parameters (such as aerosol extinction coefficient, Rayleigh scattering coefficient) through the multi-layer perceptron (MLP) structure. The output layer generates the initial atmospheric optical thickness distribution. The initial atmospheric optical thickness distribution is in terms of azimuth angle, and each azimuth angle partition corresponds to an optical thickness value, reflecting the total extinction integral of light in the vertical path in that direction. For example, if the gradient value of the 30-45 degree partition in the input optical attenuation gradient distribution is -1.5, the network may infer the presence of a high aerosol concentration in that direction and output a higher optical thickness value (such as 0.8). In the training stage, the atmospheric scattering correction network uses the true atmospheric optical thickness values (such as measured by a sun photometer) annotated in the historical dataset and the corresponding optical attenuation gradient distribution for supervised learning to optimize the weight parameters of the multi-layer perceptron to minimize the mean square error between the predicted value and the true value. Through this step, the optical attenuation gradient distribution is converted into an initial atmospheric optical thickness distribution with physical significance, providing a basis for subsequent time correction.

[0041] Step S240: Obtain the solar altitude angle parameter corresponding to the shooting time, and perform angular compensation on the initial atmospheric optical thickness distribution through the time correction module in the atmospheric scattering correction network to generate a corrected optical thickness distribution aligned with the solar azimuth.

[0042] The solar altitude angle parameter refers to the angle between the sun's rays and the horizon at the shooting moment (e.g., 90 degrees at noon and 0 degrees at sunrise), and its value is calculated through astronomical algorithms (such as the solar position model based on longitude, latitude, and time). The time correction module is an accessory structure of the atmospheric scattering correction network, which receives the initial atmospheric optical thickness distribution and the solar altitude angle parameter, and performs angular compensation on the initial distribution through linear transformation or non-linear activation functions. Specifically, when the solar altitude angle is low (such as in the early morning or evening), the sunlight penetrates a longer path through the atmosphere, resulting in a higher estimated value of the optical thickness; the time correction module dynamically adjusts the initial optical thickness value according to the solar altitude angle. For example, when the solar altitude angle is 30 degrees, the initial value is multiplied by the correction factor cos(30°)=0.866 to compensate for the influence of the path length. After angular compensation, the corrected optical thickness distribution is aligned with the actual solar azimuth, eliminating the systematic deviation caused by the change of the sun's position. For example, if the value of the initial atmospheric optical thickness distribution at a certain azimuth angle is 0.8 and the solar altitude angle is 45 degrees, the corrected value is adjusted to 0.8×cos(45°)=0.566, which more accurately reflects the atmospheric extinction characteristics in the vertical direction.

[0043] Step S250: Perform spatial projection of the corrected optical thickness distribution on the fan-shaped observation sub-regions in the target area to generate local atmospheric optical thickness parameters corresponding to each fan-shaped observation sub-region.

[0044] Spatial projection is to map the azimuth angle dimension parameters of the corrected optical thickness distribution to the fan-shaped observation sub-regions pre-divided in the target area. For example, if the corrected optical thickness distribution is divided into 24 azimuth angle partitions at 15-degree intervals, and the fan-shaped observation sub-regions in the target area are also divided at 15-degree intervals, then the optical thickness parameter of each partition directly corresponds to the corresponding value in the corrected distribution. If the division interval of the fan-shaped observation sub-regions is finer (such as 5 degrees), then the 15-degree interval parameters in the corrected distribution need to be converted to 5-degree interval parameters through linear interpolation. For example, for the 10-15 degree fan-shaped sub-region, its optical thickness parameter can be interpolated from the parameters of the 0-15 degree and 15-30 degree partitions in the corrected distribution. Finally, each fan-shaped observation sub-region is assigned a local atmospheric optical thickness parameter, forming a set of characteristic parameters that are spatially aligned with the lidar scan data. For example, the "30-45 degree fan-shaped sub-region" corresponds to the optical thickness parameter 0.566.

[0045] As an implementation manner, in step S200, path feature encoding is performed on the lidar scan data to generate extinction coefficient distribution parameters, which may specifically include:

[0046] Step S260: Segment the aerosol extinction coefficient distribution data on the horizontal scan path according to a preset distance resolution to generate a plurality of radial distance segments corresponding to the fan-shaped observation sub-regions.

[0047] The aerosol extinction coefficient distribution data on the horizontal scan path is collected by the lidar along the radial path at a preset distance resolution (such as one sampling point every 50 meters). The preset distance resolution determines the division length of the radial distance segments. For example, if the distance resolution is 50 meters and the maximum detection distance is 5000 meters, the horizontal scan path is divided into 100 radial distance segments at 50-meter intervals (0 - 50 meters, 50 - 100 meters... 4950 - 5000 meters). Each radial distance segment corresponds to a continuous spatial unit within the angular range (such as 30 - 45 degrees) of a fan-shaped observation sub-region. For example, for the 30 - 45 degree fan-shaped sub-region, the corresponding radial distance segment includes multiple intervals such as 0 - 50 meters, 50 - 100 meters within this angular range, and each interval contains the extinction coefficient sampling points of the lidar within this angular and distance range. Through this step, the lidar scan data is structurally segmented, providing input for subsequent statistical feature extraction.

[0048] Step S270: For each radial distance segment, extract the continuous extinction coefficient sampling points within the distance segment in the aerosol extinction coefficient distribution data, and calculate the mean extinction coefficient of the continuous extinction coefficient sampling points.

[0049] The continuous extinction coefficient sampling points refer to the extinction coefficient measurement values at adjacent distance positions within the same radial distance segment. For example, for the 50 - 100 meter distance segment of the 30 - 45 degree fan-shaped sub-region, the lidar may collect extinction coefficient values at 50 meters, 60 meters... 100 meters within this interval (such as 0.1 / km, 0.12 / km... 0.15 / km). The mean extinction coefficient is obtained by arithmetically averaging the extinction coefficient values of all sampling points within this distance segment. For example, the mean of the above sampling points is (0.1 + 0.12 +... + 0.15) / N, where N is the number of sampling points. The mean extinction coefficient reflects the overall level of aerosol concentration within this distance segment. A high mean indicates low visibility, while a low mean corresponds to high visibility.

[0050] Step S280: Generate the extinction integral contribution value of each radial distance segment according to the mean extinction coefficient and the length of the corresponding radial distance segment.

[0051] The extinction integral contribution value is calculated by multiplying the mean extinction coefficient by the length of the radial distance segment, which characterizes the cumulative impact of this distance segment on the overall visibility calculation. For example, if the mean extinction coefficient of a certain radial distance segment is 0.12 / km and the length is 50 meters (0.05 km), then the extinction integral contribution value is 0.12×0.05 = 0.006. This value quantifies the total attenuation of light when propagating through this distance segment. The higher the contribution value, the stronger the restrictive effect of this area on visibility. By accumulating the extinction integral contribution values of all radial distance segments, the total extinction integral along the lidar scanning path can be obtained, which is used for visibility conversion.

[0052] Step S290: Arrange the extinction integral contribution values in the spatial order of the radial distance segments to generate a sequence of extinction integral values on the horizontal scanning path.

[0053] The spatial order refers to the arrangement order of the radial distances along the lidar scanning path from near to far. For example, for the 30 - 45 degree fan-shaped sub-region, the extinction integral contribution values are arranged in the order of 0 - 50 meters, 50 - 100 meters... 4950 - 5000 meters, forming a sequence of extinction integral values containing 100 elements. This sequence reflects the trend of visibility attenuation with distance. For example, the contribution values in the short-distance segments are relatively low (low aerosol concentration), and the contribution values in the long-distance segments gradually increase (accumulation of aerosol concentration).

[0054] Step S2100: Based on the sequence of extinction integral values and the preset visibility conversion relationship, map the sequence of extinction integral values to the extinction coefficient distribution parameters that match the fan-shaped observation sub-region.

[0055] The preset visibility conversion relationship is a mathematical correlation between the extinction integral and visibility established based on an atmospheric physical model (such as the Koschmieder formula). The Koschmieder formula defines the visibility V = 3.912 / σ, where σ is the extinction coefficient. However, since the sequence of extinction integral values reflects the cumulative attenuation along the path, it is necessary to convert it to the equivalent extinction coefficient distribution parameters through an integral inversion method. For example, if the total extinction integral of a certain fan-shaped sub-region is Σ(σ i ×L i ) = 0.5 (σ i is the mean extinction coefficient of each distance segment, and L i is the length of the distance segment), and the total path length is 5 km, then the equivalent average extinction coefficient is 0.5 / 5 = 0.1 / km, corresponding to the visibility V = 3.912 / 0.1 ≈ 39.12 km. Through this conversion, the sequence of extinction integral values is mapped to the extinction coefficient distribution parameters of each fan-shaped sub-region, forming a multi-modal input that is spatially aligned with the output of image feature encoding (local atmospheric optical thickness parameters) for subsequent processing by the fusion model.

[0056] As an implementation manner, in step S300, the atmospheric optical thickness parameter is corrected in the horizontal direction by a multi-modal fusion model to generate a horizontal optical thickness correction parameter, which may specifically include:

[0057] Step S310: Input the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into the feature alignment layer of the multi-modal fusion model to generate a first feature vector of the atmospheric optical thickness parameter in the horizontal angle and a second feature vector of the extinction coefficient distribution parameter in the distance dimension.

[0058] The feature alignment layer of the multi-modal fusion model is used to convert feature parameters of different modalities into vector representations of a unified dimension to achieve spatial alignment and interaction of cross-modal features. The atmospheric optical thickness parameter takes the fan-shaped observation sub-region of the target area as the basic unit, and each unit corresponds to the local atmospheric optical thickness value within an azimuth range (such as a fan-shaped sub-region with a 15-degree interval), forming a one-dimensional parameter sequence indexed by the horizontal angle. For example, if the target area is divided into 24 fan-shaped sub-regions with a 15-degree interval, the atmospheric optical thickness parameter sequence contains 24 elements, and each element represents the vertical atmospheric extinction integral value of the corresponding fan-shaped sub-region. The feature alignment layer first unfolds this sequence into a first feature vector in the horizontal angle dimension, where the value of each dimension directly corresponds to the local atmospheric optical thickness parameter of the fan-shaped sub-region. At the same time, the extinction coefficient distribution parameter takes the radial distance segment of the lidar scanning path as the basic unit, and each unit corresponds to the equivalent extinction coefficient value within a specific distance interval (such as a radial distance segment with a 50-meter interval), forming a one-dimensional parameter sequence indexed by the distance. For example, if the maximum detection distance is 5000 meters and the distance resolution is 50 meters, the extinction coefficient distribution parameter sequence contains 100 elements, and each element represents the equivalent extinction coefficient after conversion of the extinction integral contribution value of the corresponding distance segment. The feature alignment layer unfolds this sequence into a second feature vector in the distance dimension, where the value of each dimension corresponds to the extinction coefficient distribution parameter of the radial distance segment. Through this step, two types of heterogeneous feature parameters are converted into vector forms with matching dimensions. The first feature vector takes the horizontal angle as the coordinate axis (such as dimensions 1 to 24 corresponding to fan-shaped sub-regions from 0 - 15 degrees to 345 - 360 degrees), and the second feature vector takes the distance as the coordinate axis (such as dimensions 1 to 100 corresponding to distance segments from 0 - 50 meters to 4950 - 5000 meters), providing a standardized input for subsequent cross-modal association.

[0059] Step S320: Calculate the association weights between the first feature vector and the second feature vector through the cross-attention mechanism of the multi-modal fusion model, and perform weighted fusion on the first feature vector based on the association weights to generate a fused feature vector.

[0060] The cross-attention mechanism is the core module of the multi-modal fusion model. Its function is to model the spatial correlation between features of different modalities and achieve feature complementarity through dynamic weight allocation. Specifically, the cross-attention mechanism takes the first feature vector as the query vector (Query), the second feature vector as the key vector (Key) and the value vector (Value), and calculates the correlation weight between the query vector and the key vector through matrix operations. For example, for a certain sector sub-region feature in the first feature vector (such as the atmospheric optical thickness parameter in the 60-75 degree sector sub-region corresponding to dimension 5), the cross-attention mechanism traverses all distance segment features in the second feature vector (such as the extinction coefficient distribution parameters from 0-50 meters to 4950-5000 meters corresponding to dimensions 1 to 100), calculates the similarity score between the two, and normalizes it to the correlation weight through the Softmax function. The correlation weight reflects the spatial dependence strength between the atmospheric optical thickness parameter of this sector sub-region and the extinction coefficient distribution parameters of each distance segment. For example, if the atmospheric optical thickness parameter in the 60-75 degree sector sub-region has a high correlation weight (such as 0.15) with the extinction coefficient distribution parameter in the 2000-2050 meter distance segment, it indicates that there is a significant correlation between the vertical extinction characteristics of this sector sub-region and the horizontal extinction characteristics in the 2000-2050 meter distance segment. Based on the correlation weight, a weighted sum of the second feature vector is performed to generate a context-enhanced feature for the first feature vector, and then a residual connection is made with the original first feature vector to form a fused feature vector. The fused feature vector retains the spatial distribution characteristics of the original atmospheric optical thickness parameter and at the same time incorporates the horizontal distance details provided by the lidar extinction coefficient distribution parameter, thereby enhancing the completeness of feature expression.

[0061] Step S330: Input the fused feature vector into the fully-connected layer of the multi-modal fusion model to generate a horizontal optical thickness correction parameter.

[0062] The fully connected layer is the output module of the multi-modal fusion model, which maps the high-dimensional fusion feature vector to the target parameter space through a multi-layer perceptron structure. The fusion feature vector contains the feature information of the horizontal angle dimension (such as 24 fan-shaped sub-regions) and the context information introduced by cross-modal association. The fully connected layer first compresses the dimension of the input vector to the same scale as the target parameter through a linear transformation (such as the 24-dimensional output corresponding to the correction parameters of 24 fan-shaped sub-regions), and then performs a non-linear correction on the output through an activation function (such as ReLU or Sigmoid) to ensure that the correction parameters are within a physically reasonable numerical range. For example, if the eigenvalue of a certain fan-shaped sub-region in the fusion feature vector is mapped to 0.75 after being processed by the fully connected layer, this value represents the corrected horizontal optical thickness parameter, and its value may be higher or lower than the original atmospheric optical thickness parameter (such as the original value being 0.6), reflecting the correction effect of the lidar extinction coefficient distribution parameter on the vertical extinction estimation of this fan-shaped sub-region. Finally, the fully connected layer outputs a sequence of horizontal optical thickness correction parameters, each parameter corresponding to the corrected vertical extinction integral value of a fan-shaped sub-region. Based on retaining the atmospheric vertical extinction information, this parameter sequence fuses the high-resolution extinction coefficient data of the lidar horizontal scan to achieve precise correction of the optical thickness in the horizontal direction. For example, if the estimated value of the original atmospheric optical thickness parameter of a certain fan-shaped sub-region is too high due to cloud interference, and the lidar data shows that the extinction coefficient in the short-distance section in this direction is low, then the corrected parameter will be appropriately lowered to more accurately reflect the actual atmospheric extinction characteristics.

[0063] As an implementation manner, the training process of the above multi-modal fusion model may include the following steps:

[0064] Step S10: Obtain a historical multi-modal training data set, which includes historical image observation data, historical lidar scan data, and the corresponding true visibility distribution map.

[0065] The historical multi-modal training dataset consists of multiple sets of spatio-temporally aligned observational data and ground truth labels, which are used to train the parameters of the multi-modal fusion model. The historical image observational data refers to the collection of visible light images of the target area captured by visible light imaging devices during a certain past period. Each image is attached with an accurate timestamp (such as "2022-05-10 14:30:00 UTC") and spatial position information (such as longitude and latitude coordinates, azimuth angle and elevation angle parameters of the shooting). The historical lidar scan data is the distribution data of aerosol extinction coefficients collected along the horizontal scan path by the lidar device during the same period. Each scan data records the start timestamp of the scan, the horizontal scan angle range (such as 0-360 degrees), and the sequence of extinction coefficient sampling points along the radial path. The ground truth visibility distribution map is generated by professional meteorological observation equipment (such as transmissive visibility meters or multi-point scattering sensor networks) under the same spatio-temporal conditions, and stores the measured visibility values (such as horizontal visibility values in meters) of each spatial unit in the target area in the form of a grid in the Cartesian coordinate system or polar coordinate system. For example, the historical image observational data may include a visible light image taken at 14:30 on May 10, 2022. The corresponding lidar scan data is the distribution of extinction coefficients of a 0-360 degree horizontal scan at the same moment, and the ground truth visibility distribution map is the visibility value of each 500 m × 500 m grid cell measured by the sensor network at that moment. The dataset needs to ensure time synchronization (such as the timestamp deviation is less than 1 second) and spatial coverage consistency (such as the lidar scan range completely overlaps with the geographical area of the ground truth visibility distribution map) to support the accuracy of the supervision signal for model training.

[0066] Step S20: Extract features from the historical image observational data to generate historical atmospheric optical thickness parameters, and extract features from the historical lidar scan data to generate historical extinction coefficient distribution parameters.

[0067] The feature extraction process of historical image observation data is consistent with steps S200 - S250, and historical atmospheric optical thickness parameters are generated through image feature encoding. Specifically, for each historical visible light image, first, the pixel brightness distribution sequence of its sky region is extracted (such as the brightness mean values of 24 sector partitions divided at 15 - degree intervals), the optical attenuation gradient distribution in the horizontal direction is calculated, and then it is input into a pre - trained atmospheric scattering correction network. Combining with the solar altitude angle parameter corresponding to the shooting time, the corrected optical thickness distribution is generated and projected onto the sector - shaped observation sub - regions to obtain the historical atmospheric optical thickness parameters of each sub - region. For example, for the visible light image at 14:30 on May 10, 2022, after processing, a sequence of historical atmospheric optical thickness parameters of 24 sector - shaped sub - regions is generated, and each parameter represents the vertical extinction integral value within the corresponding azimuth range (such as 30 - 45 degrees). The feature extraction process of historical lidar scan data is consistent with steps S260 - S2100. The scan path is segmented according to a preset distance resolution (such as 50 meters), the mean extinction coefficient and the extinction integral contribution value of each radial distance segment are calculated, and they are mapped to historical extinction coefficient distribution parameters through the visibility conversion relationship. For example, for the lidar scan data at the same moment, after processing, a sequence of historical extinction coefficient distribution parameters of 100 distance segments (from 0 - 50 meters to 4950 - 5000 meters) is generated, and each parameter represents the equivalent extinction coefficient value of the corresponding distance segment.

[0068] Step S30: Input the historical atmospheric optical thickness parameters and the historical extinction coefficient distribution parameters into the initial multi - modal fusion model to generate predicted visibility distribution features.

[0069] The initial multi - modal fusion model is an untrained neural network structure, and its architecture is the same as that of the model in the inference stage, including a feature alignment layer, a cross - attention mechanism, and a fully - connected layer. The historical atmospheric optical thickness parameters and the historical extinction coefficient distribution parameters are input into the model in a spatially aligned manner: the atmospheric optical thickness parameters are arranged as feature vectors in the horizontal angle dimension according to the sector - shaped sub - regions (such as 24 - dimensional vectors), and the extinction coefficient distribution parameters are arranged as feature vectors in the distance dimension according to the distance segments (such as 100 - dimensional vectors). The feature alignment layer converts the two types of vectors into an intermediate representation of a unified dimension (such as mapping to a 128 - dimensional vector through a linear transformation), the cross - attention mechanism calculates the correlation weights between the two, generates a fused feature vector, and the fully - connected layer further maps the fused feature to the predicted visibility distribution features. The predicted visibility distribution features are represented in the form of a grid in the polar coordinate system. For example, it includes 2400 grid cells of 24 angular partitions (every 15 degrees) and 100 distance segments (every 50 meters), and each cell stores the predicted visibility value. For example, inputting the atmospheric optical thickness parameter of 0.6 in the 30 - 45 - degree sector - shaped sub - region and the extinction coefficient distribution parameter of 0.12 / km in the corresponding distance segment, the model may output a predicted visibility of 8000 meters for the distance segment of 500 - 1000 meters in this region.

[0070] Step S40: Calculate the spatial consistency loss value between the predicted visibility distribution feature and the true visibility distribution map, and adjust the parameters of the initial multi-modal fusion model based on the spatial consistency loss value until the spatial consistency loss value converges.

[0071] The spatial consistency loss value is used to quantify the spatial matching degree between the predicted visibility distribution feature and the true visibility distribution map. Its calculation process includes polar coordinate transformation, difference calculation, and normalization processing. First, the predicted visibility distribution feature and the true visibility distribution map are respectively transformed into the first visibility matrix and the second visibility matrix in the polar coordinate system (Steps S41 - S415). Subsequently, the square of the visibility difference between the two is calculated for each grid cell (Step S42). Finally, all the differences are normalized and summed (Step S43). For example, if the true visibility of a certain grid cell is 10,000 meters and the predicted value is 9,500 meters, the difference is 500 meters, and the squared difference is 250,000. The sum of the squared differences for all 2,400 grid cells is divided by the total number of grids to obtain the mean squared error (MSE) as the spatial consistency loss value. The model is trained using the backpropagation algorithm. The gradient is calculated based on the loss value and the network weights are updated (such as using the Adam optimizer) until the fluctuation range of the loss value is less than a preset threshold (such as the loss value change is less than 0.001) in multiple consecutive training epochs, which is regarded as convergence.

[0072] As an implementation, Step S40, calculating the spatial consistency loss value between the predicted visibility distribution feature and the true visibility distribution map, may specifically include:

[0073] Step S41: Respectively transform the predicted visibility distribution feature and the true visibility distribution map into the first visibility matrix and the second visibility matrix in the polar coordinate system.

[0074] The polar coordinate transformation aims to unify the spatial representation forms of the prediction result and the true value, facilitating cell-by-cell comparison. The dimensions of the first visibility matrix and the second visibility matrix are defined by the polar coordinate grid parameters, including the angle dimension division sequence (such as 0 degrees, 15 degrees... 345 degrees) and the distance dimension division sequence (such as 0 meters, 50 meters... 5000 meters). Each visibility value in the predicted visibility distribution feature is converted into an angle value (arctan(y / x) = 53.13 degrees) and a distance value (√(x 2 + y 2) = 5000 m), and map it to the corresponding cell according to the grid parameters (for example, 53.13 degrees belongs to the 45 - 60 degree angle partition, and 5000 m belongs to the 4950 - 5000 m distance segment). The conversion process of the true visibility distribution map is the same, but its original coordinates may already be in polar coordinates (such as the sensor directly outputs polar coordinate data), and in this case, only resampling is required according to the same grid parameters.

[0075] As an implementation, in step S41, convert the predicted visibility distribution characteristics and the true visibility distribution map into a first visibility matrix and a second visibility matrix in the polar coordinate system respectively, including:

[0076] Step S411: Extract the spatial coverage range of the predicted visibility distribution characteristics, and generate polar coordinate grid parameters based on the maximum radial distance of the spatial coverage range and the preset angular resolution. The polar coordinate grid parameters include the angular dimension division sequence and the distance dimension division sequence.

[0077] The spatial coverage range is determined by the maximum radial distance (such as 5000 m) and the angular coverage range (such as 0 - 360 degrees) of the predicted visibility distribution characteristics. The preset angular resolution determines the division interval of the angular dimension (for example, 15 degrees corresponds to 24 partitions), and the distance dimension division sequence is generated according to the preset distance resolution (such as 50 m), arranged at equal intervals from 0 m to the maximum radial distance. For example, the maximum radial distance of 5000 m and the distance resolution of 50 m generate 100 distance segments, and the angular resolution of 15 degrees generates 24 angular partitions. The polar coordinate grid parameters are defined as the angular dimension division sequence [0°, 15°, 30°…345°] and the distance dimension division sequence [0 m, 50 m, 100 m…5000 m].

[0078] Step S412: Convert each visibility value in the predicted visibility distribution characteristics into an angular value and a distance value in the polar coordinate system according to its spatial coordinates in the Cartesian coordinate system, and map the angular value and the distance value to the corresponding polar coordinate grid cell according to the angular dimension division sequence and the distance dimension division sequence.

[0079] The conversion formula from Cartesian coordinates to polar coordinates is: angular value θ = arctan(y / x) (adjusted to the 0 - 360 degree range), distance value r = √(x 2 + y 2 ). For example, the Cartesian coordinates of a certain cell in the predicted visibility distribution characteristics are (2121 m, 2121 m), the converted angular value is 45 degrees, and the distance value is 3000 m (√(2121 2 + 2121 2) = 3000 m), mapping to the grid cell corresponding to the angular partition of 45 - 60 degrees and the distance segment of 2950 - 3000 m in the polar coordinate grid. If the angular resolution of the grid parameters is 15 degrees and the distance resolution is 50 m, this cell is classified into the cross - grid of the angular partition of 45 - 60 degrees and the distance segment of 2950 - 3000 m.

[0080] Step S413: For the blank cells in the polar coordinate grid cells that are not covered by the predicted visibility distribution characteristics, perform radial interpolation filling based on the visibility values of adjacent grid cells to generate the first visibility matrix covering the complete polar coordinate grid.

[0081] A blank cell refers to a cell in the polar coordinate grid without a mapped predicted visibility value, which may be caused by the sparsity of the original prediction data or coordinate conversion errors. Radial interpolation filling is performed along the distance dimension. For example, for a blank cell (angular partition 30 - 45 degrees, distance segment 2050 - 2100 m), find the nearest non - blank cells within the same angular partition (such as the cells of 2000 - 2050 m and 2100 - 2150 m), and calculate the linear interpolation result of their visibility values. If the visibility of the 2000 - 2050 m cell is 8000 m and that of the 2100 - 2150 m cell is 7500 m, the interpolation result of the blank cell is (8000 + 7500) / 2 = 7750 m. By this method, it is ensured that the first visibility matrix has no missing values, which is convenient for subsequent loss calculation.

[0082] Step S414: Synchronously convert each true visibility value in the true visibility distribution map into the true angular value and true distance value in the polar coordinate system according to its spatial coordinates, and map the true angular value and true distance value to the corresponding polar coordinate grid cells based on the same polar coordinate grid parameters.

[0083] The spatial coordinates of the true visibility distribution map may be in the Cartesian coordinate system (such as UTM coordinates) or the polar coordinate system, and need to be uniformly converted to the same polar coordinate grid parameters as the predicted visibility distribution characteristics. For example, the Cartesian coordinates of a cell in the true visibility distribution map are (1500 m, 1500 m), which are converted to an angular value of 45 degrees and a distance value of 2121 m in the polar coordinate system, and mapped to the grid cell corresponding to the angular partition of 45 - 60 degrees and the distance segment of 2100 - 2150 m. If the true visibility of this cell is 8500 m, the corresponding grid cell of the second visibility matrix is assigned 8500 m.

[0084] Step S415: For the blank cells in the polar coordinate grid cells that are not covered by the true visibility values, perform filling using the same interpolation rule as the radial interpolation filling to generate the second visibility matrix covering the complete polar coordinate grid.

[0085] To ensure a one-to-one correspondence between the grid cells of the first visibility matrix and the second visibility matrix, the blank cells in the true value distribution map need to be filled using the same interpolation method as the predicted data. For example, if there is no data in the 2000 - 2050 m cell of the true value distribution in the 30 - 45 degree angular partition, linear interpolation is performed based on the true value visibilities of the 1950 - 2000 m and 2050 - 2100 m cells in the same angular partition. Through this step, the second visibility matrix completely covers all polar coordinate grid cells, with exactly the same spatial structure as the first visibility matrix, thus supporting accurate loss value calculation.

[0086] Step S42: Calculate the visibility difference between the first visibility matrix and the second visibility matrix within the same angular and distance cells.

[0087] Compare the visibility values of the first visibility matrix (predicted value) and the second visibility matrix (true value) cell by cell in the grid, and calculate the absolute difference or squared difference. For example, if the predicted value of the cell in the 30 - 45 degree angular partition and the 2000 - 2050 m distance segment of the first visibility matrix is 7800 m, and the true value of the corresponding cell in the second visibility matrix is 8000 m, then the difference is 200 m and the squared difference is 40000. Traverse all grid cells (e.g., 24 angular partitions × 100 distance segments = 2400 cells), and accumulate all differences to generate the total loss value.

[0088] Step S43: Normalize the visibility difference to generate the spatial consistency loss value.

[0089] The purpose of normalization is to eliminate the dimension difference and balance the contributions of different visibility ranges. For example, using the mean squared error (MSE) as the loss function: sum the squared differences of each grid cell and divide by the total number of grids. If the total squared difference is 9600000 m 2 , and the number of grids is 2400, then the spatial consistency loss value is 9600000 / 2400 = 4000. In addition, dynamic weights (such as higher weights for low visibility areas) can be introduced or relative errors (such as the difference divided by the true value) can be used to enhance the sensitivity of the model to key regions. The normalized loss value is used as a training signal to drive the optimization of model parameters until the prediction result and the true value distribution reach high spatial consistency.

[0090] As an implementation, in step S400, based on the horizontal optical depth correction parameter and the extinction coefficient distribution parameter, perform spatial continuity fusion to generate the spatial continuous visibility distribution feature of the target area, which can specifically include:

[0091] Step S410: Divide the horizontal optical depth correction parameter into multiple angular sub - parameters according to a preset angular resolution, and divide the extinction coefficient distribution parameter into multiple distance sub - parameters according to a preset distance resolution.

[0092] The preset angular resolution is used to define the division granularity of the horizontal optical thickness correction parameter in the angular dimension. For example, when set to 1 degree, the complete horizontal scan range (0 - 360 degrees) of the target area is divided into 360 angular sub-parameters, and each angular sub-parameter corresponds to a fan-shaped sub-area with a 1-degree interval (such as 0 - 1 degree, 1 - 2 degrees... 359 - 360 degrees). The horizontal optical thickness correction parameter, as the initial input, may have an inconsistent original resolution with the preset angular resolution (such as originally divided into 24 fan-shaped sub-areas with a 15-degree interval). At this time, the optical thickness parameters of each original fan-shaped sub-area need to be expanded at 1-degree intervals through linear interpolation to generate 360 angular sub-parameters. For example, if the optical thickness correction parameter of the original 0 - 15-degree fan-shaped sub-area is 0.6, it is evenly distributed to 15 angular sub-parameters of 0 - 1 degree, 1 - 2 degrees... 14 - 15 degrees through interpolation, and the value of each sub-parameter is 0.6. Similarly, the extinction coefficient distribution parameter is divided into multiple distance sub-parameters according to the preset distance resolution (such as 10 meters). For example, when the maximum radial distance is 5000 meters, the original extinction coefficient distribution parameter is divided into 100 distance segments at 50-meter intervals, and now it is subdivided into 500 distance sub-parameters at 10-meter intervals (such as 0 - 10 meters, 10 - 20 meters... 4990 - 5000 meters) through interpolation, and each sub-parameter inherits the extinction coefficient value of the corresponding original distance segment. Through this step, both types of parameters achieve high-resolution division in the angular and distance dimensions, providing refined input for subsequent spatial continuity fusion.

[0093] Step S420: Perform linear interpolation on each angular sub-parameter and the corresponding distance sub-parameter to generate the visibility interpolation result within the continuous angular-distance grid.

[0094] The angle - distance grid is formed by the intersection of the angle dimension (0 - 360 degrees, 1 - degree interval) and the distance dimension (0 - 5000 meters, 10 - meter interval). Each grid cell corresponds to the visibility value within a specific angle range and distance range. The angular sub - parameters (360 in number) of the horizontal optical thickness correction parameter and the distance sub - parameters (500 in number) of the extinction coefficient distribution parameter are used to generate the visibility value of each grid cell through bilinear interpolation. Specifically, for a certain grid cell (such as an angle of 5 - 6 degrees and a distance of 100 - 110 meters), the visibility values of its neighboring angular sub - parameters (5 degrees and 6 degrees) and distance sub - parameters (100 meters and 110 meters) are found, and the weighted average in the angular and distance directions is calculated. For example, if the horizontal optical thickness correction parameter at 5 degrees is 0.62 and the extinction coefficient distribution parameter at 100 meters is 0.1 / km, and the horizontal optical thickness correction parameter at 6 degrees is 0.63 and the extinction coefficient distribution parameter at 110 meters is 0.09 / km, then the visibility value of this grid cell is jointly determined by the two: visibility value = horizontal optical thickness correction parameter × extinction coefficient distribution parameter, that is, 0.62×0.1 = 0.062 and 0.63×0.09 = 0.0567. Then, interpolation is performed according to the distance weight (the distance difference between 100 meters and 110 meters is 10 meters, and the proportion of the current cell's distance from 100 meters is 0.9), and the final value is 0.062×0.9 + 0.0567×0.1≈0.061. After traversing all grid cells, an interpolation result matrix containing 360×500 = 180,000 visibility values is generated, covering the discrete visibility distribution of the entire space of the target area. Although this interpolation result covers the complete grid, local mutations may occur due to sensor noise or interpolation errors and further optimization is required.

[0095] Step S430: Input the visibility interpolation result into the spatial continuity optimization network, and smooth the visibility interpolation results of adjacent grids through the spatial continuity optimization network to generate the spatial - continuous visibility distribution characteristics.

[0096] The spatial continuity optimization network is a deep learning model based on a convolutional neural network (CNN). Its input layer receives the visibility interpolation result matrix (360×500). The hidden layer extracts the spatial correlation features of adjacent grid cells through multiple sets of convolutional kernels, and the output layer generates a smoothed visibility distribution matrix. For example, the first convolutional layer uses a 3×3 convolutional kernel, slides along the angular and distance directions, calculates the weighted sum of each grid cell and its adjacent 8 cells, and outputs a feature map through an activation function (such as ReLU); the second convolutional layer further aggregates context information in a larger range, and finally the output layer maps the feature map to the visibility value space through 1×1 convolution. During the training process, the spatial continuity optimization network learns to suppress high-frequency noise in the interpolation results (such as abnormally high values in a single grid) and enhance spatial continuity (such as the gradual change of visibility in adjacent angles and distances). For example, if the visibility of the cell at angle 5 degrees and distance 100 meters in the interpolation result is 0.061, while the visibility of the adjacent cell at angle 5 degrees and distance 110 meters is 0.12 (a sudden change due to noise), the network may adjust the latter to 0.065 through convolution operations to make it smoothly transition with the surrounding cells. The finally output spatially continuous visibility distribution features eliminate local discontinuities while retaining the trend of the original data, forming a physically reasonable high-resolution visibility distribution.

[0097] Among them, the training process of the above-mentioned spatial continuity optimization network may include the following steps:

[0098] Step S4301: Obtain synthetic visibility training data, which includes simulated visibility interpolation results and corresponding continuous visibility distribution labels.

[0099] The synthetic visibility training data is generated by numerical simulation or adding noise to real data. The numerical simulation method generates idealized continuous visibility distribution labels based on an atmospheric physics model (such as the Monte Carlo radiative transfer model), and then artificially adds random noise (such as Gaussian noise or impulse noise) to generate simulated visibility interpolation results. For example, the visibility in a certain area of the continuous visibility distribution label decreases linearly from 5000 meters (angles 0 - 90 degrees) to 1000 meters (angles 270 - 360 degrees), and the simulated interpolation result superimposes noise with a standard deviation of 500 meters on this label to form a training input containing abnormal fluctuations and discontinuities. The real data noise addition method uses the historical true visibility distribution map as a label, downsamples it (such as from a 10-meter resolution to 50 meters) and then interpolates it back to a high resolution to simulate low-quality interpolation results. The synthetic data needs to cover various atmospheric conditions (such as fog, haze, clear sky) and noise patterns to ensure the generalization ability of the model.

[0100] Step S4302: Input the simulated visibility interpolation results into the initial spatial continuity optimization network to generate a predicted continuous visibility distribution.

[0101] The initial spatial continuity optimization network is an untrained CNN model, and its structure is the same as that in the inference stage. The simulated visibility interpolation results are input into the network in a grid format (such as a 360×500 matrix). Through the convolutional layers, features are extracted step by step, and the predicted continuous visibility distribution is output. For example, when inputting an interpolated result with noise (the visibility of the unit at an angle of 5 degrees and a distance of 100 meters is 0.12), the network may output a corrected visibility of 0.065 by weighting the values of the surrounding units (angles 4 - 6 degrees, distances 90 - 110 meters) with convolutional kernels. In the initial stage of training, the prediction results may deviate significantly from the labels, but they are gradually optimized under the guidance of the loss function.

[0102] Step S4303: Calculate the gradient difference loss value between the predicted continuous visibility distribution and the label of the continuous visibility distribution, and adjust the parameters of the initial spatial continuity optimization network based on the gradient difference loss value.

[0103] The gradient difference loss value is calculated by comparing the first-order gradients of the prediction results and the labels in the angular and distance directions (i.e., the visibility change rate of adjacent units). Specifically, calculate the mean square error of the horizontal angular gradient (the visibility difference between adjacent angular units) and the radial distance gradient (the visibility difference between adjacent distance units) of the predicted distribution and the label distribution respectively. For example, if the visibility difference between the unit at an angle of 5 degrees and a distance of 100 meters and the unit at an angle of 6 degrees and a distance of 100 meters in the label is -0.002, and the corresponding difference in the predicted distribution is 0.005, then the angular gradient error is (0.005 - (-0.002)) 2 = 0.000049; Accumulate all the angular and distance gradient errors and normalize them to obtain the gradient difference loss value. Through the backpropagation algorithm, the loss value drives the network to adjust the convolutional kernel weights so that the predicted distribution output by it is consistent with the label in terms of gradient distribution, thereby enforcing spatial continuity. The training termination condition is that the gradient difference loss value does not decrease significantly for multiple consecutive epochs on the validation set (such as the change amplitude is less than 0.1%), indicating that the model has fully learned the smoothing rules.

[0104] As an implementation, in step S500, converting the spatially continuous visibility distribution feature into the visibility prediction result in the polar coordinate system may specifically include:

[0105] Step S510: Generate polar coordinate grid division parameters based on the maximum radial distance of the spatially continuous visibility distribution feature and the preset angular interval of the fan-shaped observation sub-region. The polar coordinate grid division parameters include the angular axis division sequence and the distance axis division sequence.

[0106] The maximum radial distance of the spatial continuous visibility distribution characteristics is determined by the maximum effective detection range of lidar scanning, such as 5000 meters, indicating the upper limit of the radial distance from the center to the periphery of the target area covered by the visibility prediction result. The preset angular interval of the fan-shaped observation sub-region defines the division accuracy of the polar coordinate grid in the angular dimension. For example, a 15-degree interval divides the horizontal circle into 24 fan-shaped sub-regions (0-15 degrees, 15-30 degrees... 345-360 degrees). Based on the above parameters, the angular axis division sequence in the polar coordinate grid division parameters is generated by the preset angular interval, that is, [0°, 15°, 30°... 345°], and the distance axis division sequence is equally spaced from 0 meters to the maximum radial distance according to the preset distance resolution (such as 50 meters), generating a sequence of [0m, 50m, 100m... 5000m]. For example, when the maximum radial distance is 5000 meters and the distance resolution is 50 meters, the distance axis division sequence contains 100 distance segments, and each distance segment corresponds to a radial interval of 50 meters. The angular axis and the distance axis division sequences together constitute the index framework of the polar coordinate grid, ensuring the standardization of the spatial structure of the visibility prediction result.

[0107] Step S520: Extract the original coordinates of each visibility value in the spatial continuous visibility distribution characteristics in the Cartesian coordinate system, and convert the original coordinates into the angular component and the distance component in the polar coordinate system.

[0108] The visibility values of the spatial continuous visibility distribution characteristics are stored as two-dimensional grid data in the Cartesian coordinate system, and each visibility value is associated with a coordinate point (such as x = 3000 meters, y = 4000 meters). Calculate the angular component θ and the distance component r in the polar coordinate system through the coordinate conversion formula: the angular component θ = arctan(y / x) × (180 / π), and adjust it to the range of 0-360 degrees according to the coordinate quadrant; the distance component r = √(x 2 + y 2 ). For example, the Cartesian coordinates (3000 meters, 4000 meters) are converted into the angular component θ = arctan(4000 / 3000) ≈ 53.13 degrees and the distance component r = √(3000 2 + 4000 2 ) = 5000 meters. During the conversion process, it is necessary to ensure the continuity of the angular component. For example, when x is negative, the angular component needs to be increased by 180 degrees to correctly reflect the actual orientation. Through this step, all the Cartesian coordinates of the visibility values are uniformly mapped into the angular and distance parameters in the polar coordinate system, forming input data that matches the polar coordinate grid division parameters.

[0109] Step S530: Map the angular component to the corresponding angular axis grid cell according to the angular axis division sequence, and map the distance component to the corresponding distance axis grid cell according to the distance axis division sequence, generating an initial polar coordinate visibility distribution mapping table.

[0110] The angular axis division sequence and the distance axis division sequence divide the polar coordinate system into discrete grid cells, and each cell is uniquely identified by its angular range (such as 30 - 45 degrees) and distance range (such as 1000 - 1050 meters). Classify the converted angular component and distance component after step S520 into the corresponding grid cells respectively: the angular component θ = 53.13 degrees belongs to the angular axis grid cell of 45 - 60 degrees, and the distance component r = 5000 meters belongs to the distance axis grid cell of 4950 - 5000 meters. After traversing the polar coordinate parameters of all visibility values, an initial polar coordinate visibility distribution mapping table is generated, where each grid cell records the set of visibility values falling into that cell. For example, if multiple visibility values are mapped to the same grid cell (such as the cell of 45 - 60 degrees, 4950 - 5000 meters), the average value or the maximum value is taken as the initial visibility value of that cell; if there is no data mapping for a certain cell, it is marked as a blank cell. There may be data sparse areas in the initial mapping table, which need to be filled by interpolation.

[0111] Step S540: Traverse each grid cell in the initial polar coordinate visibility distribution mapping table. If there is no mapped visibility value in the current grid cell, extract the visibility values in the adjacent angular direction and distance direction of the current grid cell for two-way neighboring interpolation to generate the interpolated visibility value of the current grid cell.

[0112] Two-way neighboring interpolation is used to solve the data missing problem in the initial mapping table and ensure the integrity of the polar coordinate grid. For a blank grid cell (such as the angular range of 30 - 45 degrees and the distance range of 2000 - 2050 meters), extract the visibility values of the adjacent angular cells (such as 15 - 30 degrees and 45 - 60 degrees) along the angular axis direction, and extract the visibility values of the adjacent distance cells (such as 1950 - 2000 meters and 2050 - 2100 meters) along the distance axis direction. When interpolating, the linear weighted average method is adopted: if the visibility values of the adjacent angular cells are 8000 meters (15 - 30 degrees) and 8200 meters (45 - 60 degrees), and the visibility values of the adjacent distance cells are 7900 meters (1950 - 2000 meters) and 8100 meters (2050 - 2100 meters), then the interpolation result of the blank cell is the average value of the angular direction average value (8000 + 8200) / 2 = 8100 meters and the distance direction average value (7900 + 8100) / 2 = 8000 meters, that is, 8050 meters. If there is no valid data in a certain direction, interpolation is only performed in the other direction. Through this step, all blank cells are filled, generating visibility distribution data that seamlessly covers the polar coordinate grid.

[0113] Step S550: Rearrange the visibility values of all grid cells after interpolation according to the angular axis division sequence and the distance axis division sequence to generate a coordinate distribution map of the visibility prediction result covering the complete polar coordinate range.

[0114] The coordinate distribution map of the polar coordinate visibility prediction results has the angle axis as the circumferential direction and the distance axis as the radial direction. The visibility values of the grid cells are arranged by dividing the angle axis into sequences (0°, 15°…345°) and the distance axis into sequences (0m, 50m…5000m). For example, the first column of grid cells corresponding to 0° on the angle axis stores the visibility values for distance segments from 0 - 50 meters, 50 - 100 meters…4950 - 5000 meters in sequence. The second column corresponding to 15° on the angle axis stores the visibility values for the same distance segments, until the last column at 345° on the angle axis. The coordinate distribution map can be rendered as a polar coordinate heat map or a contour map through a visualization tool, where the color gradient or contour interval represents the level of visibility values (e.g., red indicates low visibility and green indicates high visibility). The finally output coordinate distribution map can be directly used in a meteorological warning system or a traffic management platform to provide the visibility prediction results with continuous spatial distribution in the target area. For example, it can identify a fog area with visibility less than 1000 meters within the range of 2000 - 3000 meters in the northeast direction (45 - 60 degrees), while the visibility in the southwest direction (225 - 240 degrees) is generally higher than 5000 meters.

[0115] In a possible implementation, after outputting the coordinate distribution map of the visibility prediction results, it may also include the process of generating dynamic visibility warnings, which may specifically include the following steps:

[0116] Step S600: Extract the real-time visibility value of each polar coordinate grid cell in the coordinate distribution map, continuously obtain the historical visibility prediction sequence at a preset time interval, and generate a dynamic visibility change trend map.

[0117] Each polar coordinate grid cell in the coordinate distribution map stores the visibility prediction value at the current moment. For example, the real-time visibility of the cell with an angle axis of 30 - 45 degrees and a distance axis of 2000 - 2050 meters is 8000 meters. The preset time interval is used to define the acquisition frequency of historical visibility prediction data. For example, it is collected once every 5 minutes, and 12 sets of historical prediction data within the past 1 hour are continuously obtained. The dynamic visibility change trend map stores the visibility value sequences of each grid cell at different time points in a four-dimensional time - space structure (time, angle, distance, visibility value). For example, the visibility of the cell with an angle of 30 - 45 degrees and a distance of 2000 - 2050 meters is 8000 meters at time point T1 (14:00), 7800 meters at T2 (14:05)…7500 meters at T12 (15:00), forming the visibility change sequence of this cell. By stacking the data of all grid cells along the time axis, a historical visibility prediction cube covering the entire space of the target area is generated for analyzing the spatio - temporal evolution law of visibility.

[0118] Step S700: Calculate the difference between the visibility values at adjacent time intervals in the dynamic visibility change trend graph to generate a visibility fluctuation amplitude sequence for each grid cell.

[0119] The difference calculation refers to performing a subtraction operation on the visibility values at adjacent time points for the same grid cell to obtain the visibility fluctuation amplitude. For example, the visibility difference between time point T2 (14:05) and T1 (14:00) is 7800 m - 8000 m = -200 m, and the difference between time point T3 (14:10) and T2 is 7600 m - 7800 m = -200 m, forming a fluctuation amplitude sequence [-200, -200…] for this cell. The fluctuation amplitude sequence reflects the change rate and direction of visibility. A positive fluctuation amplitude indicates an increase in visibility (such as fog dissipation), and a negative value indicates a decrease in visibility (such as fog formation). The fluctuation amplitude sequences of all grid cells are arranged in chronological order to form a dynamic visibility fluctuation matrix, whose dimension is (number of time intervals - 1) × number of angular axis units × number of distance axis units.

[0120] Step S800: Based on the visibility fluctuation amplitude sequence and the preset dynamic threshold interval, mark the grid cells that exceed the dynamic threshold interval as abnormally fluctuating cells.

[0121] The preset dynamic threshold interval is used to define the reasonable range of visibility fluctuations. For example, the threshold lower limit is set to -500 m / 5 minutes and the upper limit is set to +500 m / 5 minutes. Fluctuations outside this range are considered abnormal. The marking of abnormally fluctuating cells is achieved by traversing the fluctuation amplitude sequence of each grid cell: if the fluctuation amplitude at a certain time point is lower than the threshold lower limit or higher than the upper limit, then mark this time - space cell as abnormal. For example, the fluctuation amplitude of the cell in the range of angle 30 - 45 degrees and distance 2000 - 2050 m at time point T2 is -600 m (lower than -500 m), then it is marked as abnormal. The threshold interval can be dynamically adjusted according to historical statistical results or meteorological experience. For example, a more stringent threshold (such as ±300 m) is set during the fog season.

[0122] Step S900: Cluster the abnormally fluctuating cells into abnormally fluctuating regions according to spatial continuity, and overlay the spatial boundaries of the abnormally fluctuating regions onto the coordinate distribution map to generate a dynamic early warning annotation map.

[0123] Spatial continuity clustering uses a region growing algorithm to merge adjacent abnormal fluctuation units to form continuous regions. The definition of adjacency includes proximity in terms of both angle and distance dimensions: at the same time point, if the angle difference between two abnormal units is less than a preset angle tolerance (such as 15 degrees) and the distance difference is less than a preset distance tolerance (such as 100 meters), they are considered adjacent. For example, at time point T2, the unit with an angle of 30 - 45 degrees and a distance of 2000 - 2050 meters and the unit with an angle of 45 - 60 degrees and a distance of 2050 - 2100 meters are both marked as abnormal, and with an angle difference of 15 degrees and a distance difference of 50 meters, meeting the proximity condition, they are clustered into the same abnormal fluctuation region. The region boundary is generated by extracting the outer contour of the clustered region, for example, using the convex hull algorithm to determine the minimum bounding polygon. Based on the original coordinate distribution map, the dynamic early warning annotation map marks the spatial position and time range of the abnormal fluctuation region with highlighted color blocks or contour lines. For example, a red semi - transparent overlay is used to identify that there is an abnormal low visibility fluctuation in the region with an angle of 30 - 90 degrees and a distance of 1500 - 3000 meters during the period from 14:05 to 14:20.

[0124] In a possible implementation manner, after generating the dynamic early warning annotation map, the embodiment of the present invention may further include a process of abnormal region verification and correction, which may specifically include the following steps:

[0125] Step S1000: Obtain the historical concurrent visibility data corresponding to the abnormal fluctuation region, and extract the visibility mean sequence at the same spatial position in the historical concurrent period.

[0126] The historical concurrent visibility data refers to meteorological observation data in the same period (such as the same month, time period) in the past several years. For example, if the current date is 14:00 - 15:00 on November 15, 2023, the historical concurrent data is the measured visibility values at 14:00 - 15:00 on November 15 of each year from 2020 to 2022. For each grid unit (such as an angle of 30 - 45 degrees and a distance of 2000 - 2050 meters) within the abnormal fluctuation region, extract the visibility values at the same spatial position in the historical concurrent period, and calculate its mean sequence. For example, the visibility of this unit in 2020 is 8500 meters, in 2021 is 8200 meters, and in 2022 is 8000 meters, then the historical concurrent mean is (8500 + 8200 + 8000) / 3 ≈ 8233 meters. The mean sequence reflects the typical visibility level at this position in the historical concurrent period, providing a benchmark for deviation analysis.

[0127] Step S1100: Compare the real - time visibility values with the visibility mean sequence grid - by - grid to generate a visibility deviation distribution map.

[0128] The calculation formula for the deviation degree is: Deviation degree = (Real-time visibility value - Historical average value in the same period) / Historical average value in the same period × 100%. For example, if the real-time visibility is 7500 meters and the historical average value in the same period is 8233 meters, then the deviation degree is (7500 - 8233) / 8233 × 100% ≈ -8.9%. A negative deviation degree indicates that the real-time visibility is lower than the historical level in the same period, and a positive deviation degree indicates the opposite. The deviation degree distribution map stores the deviation degree values of each unit in the form of a polar coordinate grid. For example, the deviation degree of the unit with an angle of 30 - 45 degrees and a distance of 2000 - 2050 meters is -8.9%, and the deviation degree of the unit with an angle of 45 - 60 degrees and a distance of 2050 - 2100 meters is -12.5%. The distribution map visualizes the deviation degree through color coding. For example, dark red indicates a deviation degree < -15%, light red indicates -15% to -5%, and green indicates the normal range (-5% to +5%).

[0129] Step S1200: According to the deviation direction and amplitude of each grid unit in the visibility deviation degree distribution map, mark the grid units with a deviation degree exceeding the preset threshold as persistent abnormal units.

[0130] The preset threshold includes a deviation degree amplitude threshold and a duration threshold. For example, the amplitude threshold is set to -10% (that is, the real-time visibility is more than 10% lower than the historical average value in the same period), and the duration threshold is 3 consecutive time intervals (15 minutes). The marking of persistent abnormal units needs to meet both conditions: 1) The deviation degree is lower than -10%; 2) It meets the condition at 3 or more consecutive time points. For example, the deviation degrees of the unit with an angle of 30 - 45 degrees and a distance of 2000 - 2050 meters at time points T1 - T3 are -8.9%, -12.1%, and -15.3% respectively. It is marked as a persistent abnormal only when it is lower than -10% for 3 consecutive times. This dual-threshold mechanism can filter out short-term fluctuations and focus on long-term anomalies.

[0131] Step S1300: Based on the spatial distribution density of persistent abnormal units, divide the target area into a high-confidence abnormal area and a low-confidence abnormal area, and integrate the division results into the dynamic early warning annotation map to generate a graded early warning distribution map.

[0132] The spatial distribution density calculates the proportion of the number of persistent anomaly units in each preset sub-region (such as every 30 degrees × 500 meters). For example, the sub-region with an angle of 30 - 60 degrees and a distance of 2000 - 2500 meters contains 15 grid units, among which 12 are persistent anomaly units, and the density is 80%. The high-confidence anomaly area is defined as the sub-region with a density exceeding 50%, and the low-confidence area is the area with a density of 20% - 50%. The graded early warning distribution map, based on the dynamic early warning annotation map, differentiates high / low confidence areas with different identifiers: the high-confidence area is filled with a solid red frame and dark red, and the low-confidence area is filled with a dashed orange frame and light orange. For example, the area with an angle of 30 - 60 degrees and a distance of 2000 - 2500 meters is marked as a high-confidence anomaly area, and the area with an angle of 60 - 90 degrees and a distance of 2500 - 3000 meters is marked as a low-confidence anomaly area.

[0133] In a possible implementation manner, after generating the graded early warning distribution map, the embodiment of the present invention may further include the process of real-time visibility update and feedback, which may specifically include the following steps:

[0134] Step S1400: Extract the real-time visibility value and the corresponding timestamp of the high-confidence anomaly area, and generate a visibility attenuation rate curve.

[0135] The visibility attenuation rate refers to the amount of decrease in visibility per unit time, and the calculation formula is: rate = (current visibility - previous visibility) / time interval. For example, in the high-confidence anomaly area, the visibility at time point T1 (14:00) is 8000 meters, and at T2 (14:05) it is 7500 meters, with a time interval of 5 minutes, then the attenuation rate is (7500 - 8000) / 5 = -100 meters / minute. The attenuation rate curve is a line graph with the timestamp as the horizontal axis and the rate as the vertical axis. For example, the rate from T1 - T2 is -100 meters / minute, from T2 - T3 is -150 meters / minute, and from T3 - T4 is -200 meters / minute, and the curve shows an accelerating downward trend.

[0136] Step S1500: Input the attenuation rate curve into a pre-trained visibility trend prediction model to generate a visibility prediction attenuation path within a future time window.

[0137] The visibility trend prediction model is a time series prediction model based on LSTM (Long Short-Term Memory Network). Its input is the attenuation rate sequence of the past N time points (such as the rate values of the past 6 five-minute intervals), and the output is the predicted rate of the next M time points (such as 12 intervals in the next hour). For example, for the input sequence [-100, -150, -200, -250, -300, -350] m / min, the model predicts that the rate will continue to decline to -500 m / min at the next 12 time points. The attenuation path generates future visibility values by integrating the predicted rate: for example, if the current visibility is 7500 m and the predicted rate for the next time point is -400 m / min, then the visibility after 5 minutes is 7500 - 400×5 = 5500 m, and so on to generate the visibility attenuation path for the next 1 hour.

[0138] Step S1600: Based on the spatial coverage range of the visibility prediction attenuation path, mark the boundary of the predicted attenuation area in the hierarchical early warning distribution map.

[0139] The boundary of the predicted attenuation area is generated by connecting the grid cells where the visibility at each future time point is lower than a preset threshold (such as 3000 m). For example, the predicted path shows that the visibility in the area with an angle of 30 - 60 degrees and a distance of 2000 - 2500 m will drop to 2000 m at 15:00, and will expand to an angle of 30 - 90 degrees and a distance of 1500 - 3000 m at 15:30. The boundary marking uses a gradient color band or an arrow to indicate the attenuation direction. For example, the predicted boundary at 15:00 is marked with a purple dashed line, and the predicted boundary at 15:30 is marked with a blue dashed line in the hierarchical early warning distribution map.

[0140] Step S1700: Dynamically adjust the early warning level of the dynamic early warning marking map according to the spatial overlap ratio between the boundary of the predicted attenuation area and the real-time visibility value.

[0141] The spatial overlap ratio calculates the proportion of the intersection area of the predicted area and the current abnormal area. For example, if the predicted area covers an angle of 30 - 90 degrees and a distance of 1500 - 3000 m, and the current high-confidence abnormal area is an angle of 30 - 60 degrees and a distance of 2000 - 2500 m, the overlap ratio is (30 - 60 degrees intersection 30 - 90 degrees) × (2000 - 2500 m intersection 1500 - 3000 m) / predicted area = 100%. The early warning level is upgraded according to the overlap ratio: if the ratio ≥ 70%, the early warning level is upgraded from yellow to orange; if ≥ 90% and the predicted visibility is lower than the critical value (such as 1000 m), it is upgraded to a red early warning. The adjusted dynamic early warning marking map is dynamically refreshed by colors and icons. For example, the red early warning area shows a flashing icon and a text prompt to remind the traffic control department to initiate an emergency response.

[0142] According to another aspect of the present invention, there is also provided a spatial continuous visibility prediction device 900 based on multi-modal data analysis. Please refer to Figure 3 , which includes the following component structures:

[0143] A data acquisition module 910, configured to acquire multi-modal observation data of a target area. The multi-modal observation data includes image observation data and lidar scan data. The image observation data includes visible light images of the target area and corresponding shooting time and spatial position information, and the lidar scan data includes aerosol extinction coefficient distribution data on the horizontal scan path of the target area;

[0144] A feature encoding module 920, configured to perform image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and perform path feature encoding on the lidar scan data to generate extinction coefficient distribution parameters;

[0145] A parameter correction module 930, configured to input the atmospheric optical thickness parameters and the extinction coefficient distribution parameters into a multi-modal fusion model, and perform horizontal optical thickness correction on the atmospheric optical thickness parameters through the multi-modal fusion model to generate horizontal optical thickness correction parameters;

[0146] A parameter fusion module 940, configured to perform spatial continuity fusion based on the horizontal optical thickness correction parameters and the extinction coefficient distribution parameters to generate spatial continuous visibility distribution characteristics of the target area;

[0147] A visibility prediction module 950, configured to convert the spatial continuous visibility distribution characteristics into a visibility prediction result in the polar coordinate system and output a coordinate distribution map of the visibility prediction result.

[0148] The description of the above device embodiments is similar to the description of the above method embodiments and has similar beneficial effects to the method embodiments. In some embodiments, the functions or modules included in the device provided by the present invention can be used to execute the methods described in the above method embodiments. For the technical details not disclosed in the device embodiments of the present invention, please refer to the description of the method embodiments of the present invention for understanding.

Claims

1. A method for predicting spatial continuous visibility based on multi-modal data analysis, characterized in that, The method includes: Obtaining multi-modal observation data of a target area, where the multi-modal observation data includes image observation data and lidar scan data. The image observation data includes visible light images of the target area and corresponding shooting time and spatial position information, and the lidar scan data includes aerosol extinction coefficient distribution data on the horizontal scan path of the target area; Performing image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and performing path feature encoding on the lidar scan data to generate extinction coefficient distribution parameters; Inputting the atmospheric optical thickness parameters and the extinction coefficient distribution parameters into a multi-modal fusion model, and correcting the optical thickness in the horizontal direction of the atmospheric optical thickness parameters through the multi-modal fusion model to generate horizontal optical thickness correction parameters; Performing spatial continuity fusion based on the horizontal optical thickness correction parameters and the extinction coefficient distribution parameters to generate a spatially continuous visibility distribution feature of the target area; Converting the spatially continuous visibility distribution feature into a visibility prediction result in the polar coordinate system and outputting a coordinate distribution map of the visibility prediction result.

2. The method according to claim 1, wherein The obtaining of the multi-modal observation data of the target area includes: Extracting the shooting timestamp of the image observation data and the scan timestamp of the lidar scan data, and marking the image observation data and the lidar scan data with a timestamp difference less than a preset synchronization threshold as an associated data set with time alignment; Based on the horizontal scan angle range of the lidar scan data in the associated data set, dividing the target area into multiple fan-shaped observation sub-areas at a preset angle interval, and assigning a corresponding radial distance interval to each fan-shaped observation sub-area; According to the spatial position information of the image observation data in the associated data set, extracting the pixel region in the visible light image that overlaps with the angle range of the fan-shaped observation sub-area as a local observation image block, and spatially binding the local observation image block to the corresponding radial distance interval; Performing distance interval filtering on the lidar scan data in each fan-shaped observation sub-area, retaining the aerosol extinction coefficient distribution data within the radial distance interval, and generating an extinction coefficient subset corresponding one-to-one to the local observation image block; Arranging the local observation image blocks and the extinction coefficient subsets in the angular order of the fan-shaped observation sub-areas to generate a spatially aligned multi-modal observation data sequence.

3. The method according to claim 1, characterized in that, The performing of image feature encoding on the image observation data to generate atmospheric optical thickness parameters includes: Extracting the pixel brightness distribution sequence of the sky area in the visible light image, where the pixel brightness distribution sequence includes the brightness means of the sky area at different azimuth angles; Generating the optical attenuation gradient distribution of the visible light image in the horizontal direction according to the change trend of the brightness means of adjacent azimuth angles in the pixel brightness distribution sequence; Input the optical attenuation gradient distribution into a pre-trained atmospheric scattering correction network, and map the non-linear relationship between the optical attenuation gradient distribution and the atmospheric scattering parameters through a multi-layer perception structure in the atmospheric scattering correction network to generate an initial atmospheric optical thickness distribution; Obtain the solar altitude angle parameter corresponding to the shooting time, and perform angular compensation on the initial atmospheric optical thickness distribution through a time correction module in the atmospheric scattering correction network to generate a corrected optical thickness distribution aligned with the solar azimuth; Perform spatial projection of the corrected optical thickness distribution on the fan-shaped observation sub-regions in the target area to generate local atmospheric optical thickness parameters corresponding to each fan-shaped observation sub-region; Encoding the path features of the lidar scan data to generate an extinction coefficient distribution parameter, including: Segment the aerosol extinction coefficient distribution data on the horizontal scan path according to a preset distance resolution to generate a plurality of radial distance segments corresponding to the fan-shaped observation sub-regions; For each radial distance segment, extract the continuous extinction coefficient sampling points within the distance segment in the aerosol extinction coefficient distribution data, and calculate the mean extinction coefficient of the continuous extinction coefficient sampling points; Generate the extinction integral contribution value of each radial distance segment according to the mean extinction coefficient and the length of the corresponding radial distance segment; Arrange the extinction integral contribution values in the spatial order of the radial distance segments to generate an extinction integral value sequence on the horizontal scan path; Based on the extinction integral value sequence and a preset visibility conversion relationship, map the extinction integral value sequence to an extinction coefficient distribution parameter matching the fan-shaped observation sub-region.

4. The method according to claim 1, wherein The optical thickness correction in the horizontal direction of the atmospheric optical thickness parameter by the multi-modal fusion model to generate a horizontal optical thickness correction parameter, including: Input the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into the feature alignment layer of the multi-modal fusion model to generate a first feature vector of the atmospheric optical thickness parameter in the horizontal angle and a second feature vector of the extinction coefficient distribution parameter in the distance dimension; Calculate the correlation weight between the first feature vector and the second feature vector through the cross-attention mechanism of the multi-modal fusion model, and perform weighted fusion on the first feature vector based on the correlation weight to generate a fused feature vector; Input the fused feature vector into the fully connected layer of the multi-modal fusion model to generate the horizontal optical thickness correction parameter.

5. The method according to claim 4, characterized in that, The training process of the multi-modal fusion model includes: Obtain a historical multi-modal training data set, which includes historical image observation data, historical lidar scan data, and the corresponding true visibility distribution map; Extract features from the historical image observation data to generate historical atmospheric optical thickness parameters, and extract features from the historical lidar scan data to generate historical extinction coefficient distribution parameters; Input the historical atmospheric optical thickness parameter and the historical extinction coefficient distribution parameter into an initial multi-modal fusion model to generate a predicted visibility distribution feature; Calculate the spatial consistency loss value between the predicted visibility distribution characteristics and the true visibility distribution map, and adjust the parameters of the initial multi-modal fusion model based on the spatial consistency loss value until the spatial consistency loss value converges.

6. The method according to claim 5, characterized in that, The calculation of the spatial consistency loss value between the predicted visibility distribution characteristics and the true visibility distribution map includes: Convert the predicted visibility distribution characteristics and the true visibility distribution map into a first visibility matrix and a second visibility matrix in the polar coordinate system respectively; Calculate the visibility difference between the first visibility matrix and the second visibility matrix within the same angular and distance units; Normalize the visibility difference to generate the spatial consistency loss value.

7. The method according to claim 6, wherein The conversion of the predicted visibility distribution characteristics and the true visibility distribution map into a first visibility matrix and a second visibility matrix in the polar coordinate system respectively includes: Extract the spatial coverage range of the predicted visibility distribution characteristics, and generate polar coordinate grid parameters based on the maximum radial distance of the spatial coverage range and a preset angular resolution. The polar coordinate grid parameters include an angular dimension division sequence and a distance dimension division sequence; Convert each visibility value in the predicted visibility distribution characteristics into an angular value and a distance value in the polar coordinate system according to its spatial coordinates in the Cartesian coordinate system, and map the angular value and the distance value to the corresponding polar coordinate grid cells according to the angular dimension division sequence and the distance dimension division sequence; For the blank cells in the polar coordinate grid cells that are not covered by the predicted visibility distribution characteristics, perform radial interpolation filling based on the visibility values of adjacent grid cells to generate a first visibility matrix covering the complete polar coordinate grid; Synchronously convert each true visibility value in the true visibility distribution map into a true angular value and a true distance value in the polar coordinate system according to its spatial coordinates, and map the true angular value and the true distance value to the corresponding polar coordinate grid cells based on the same polar coordinate grid parameters; For the blank cells in the polar coordinate grid cells that are not covered by the true visibility values, perform filling using the same interpolation rule as the radial interpolation filling to generate a second visibility matrix covering the complete polar coordinate grid.

8. The method according to claim 1, wherein The generation of the spatially continuous visibility distribution characteristics of the target area by performing spatial continuity fusion based on the horizontal optical depth correction parameter and the extinction coefficient distribution parameter includes: Divide the horizontal optical depth correction parameter into multiple angular sub-parameters according to a preset angular resolution, and divide the extinction coefficient distribution parameter into multiple distance sub-parameters according to a preset distance resolution; Perform linear interpolation on each angular sub-parameter and the corresponding distance sub-parameter to generate visibility interpolation results within a continuous angular-distance grid; Input the visibility interpolation results into a spatial continuity optimization network, and smooth the visibility interpolation results of adjacent grid cells through the spatial continuity optimization network to generate the spatially continuous visibility distribution characteristics; Among them, the training method of the spatial continuity optimization network includes: Obtain synthetic visibility training data, where the synthetic visibility training data includes simulated visibility interpolation results and corresponding continuous visibility distribution labels; Input the simulated visibility interpolation results into the initial spatial continuity optimization network to generate a predicted continuous visibility distribution; Calculate the gradient difference loss value between the predicted continuous visibility distribution and the continuous visibility distribution label, and adjust the parameters of the initial spatial continuity optimization network based on the gradient difference loss value.

9. The method according to claim 2 or 3, characterized in that The conversion of the spatially continuous visibility distribution feature into a visibility prediction result in the polar coordinate system includes: Generate polar coordinate grid division parameters based on the maximum radial distance of the spatially continuous visibility distribution feature and the preset angular interval of the fan-shaped observation sub-region. The polar coordinate grid division parameters include an angular axis division sequence and a distance axis division sequence; Extract the original coordinates in the Cartesian coordinate system of each visibility value in the spatially continuous visibility distribution feature, and convert the original coordinates into angular components and distance components in the polar coordinate system; Map the angular components to the corresponding angular axis grid cells according to the angular axis division sequence, and map the distance components to the corresponding distance axis grid cells according to the distance axis division sequence to generate an initial polar coordinate visibility distribution mapping table; Traverse each grid cell in the initial polar coordinate visibility distribution mapping table. If there is no mapped visibility value in the current grid cell, extract the visibility values in the adjacent angular and distance directions of the current grid cell for two-way neighboring interpolation to generate the interpolated visibility value of the current grid cell; Rearrange the visibility values of all interpolated grid cells according to the angular axis division sequence and the distance axis division sequence to generate a coordinate distribution map of the visibility prediction result covering the complete polar coordinate range.

10. A spatial continuous visibility prediction device based on multi-modal data analysis, characterized in that It includes: A data acquisition module for acquiring multi-modal observation data of the target area. The multi-modal observation data includes image observation data and lidar scan data. The image observation data includes visible light images of the target area and corresponding shooting time and spatial position information, and the lidar scan data includes aerosol extinction coefficient distribution data on the horizontal scan path of the target area; A feature encoding module for performing image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and performing path feature encoding on the lidar scan data to generate extinction coefficient distribution parameters; A parameter correction module for inputting the atmospheric optical thickness parameters and the extinction coefficient distribution parameters into a multi-modal fusion model, and performing horizontal optical thickness correction on the atmospheric optical thickness parameters through the multi-modal fusion model to generate horizontal optical thickness correction parameters; A parameter fusion module for performing spatial continuity fusion based on the horizontal optical thickness correction parameters and the extinction coefficient distribution parameters to generate a spatially continuous visibility distribution feature of the target area; A visibility prediction module for converting the spatially continuous visibility distribution feature into a visibility prediction result in the polar coordinate system and outputting a coordinate distribution map of the visibility prediction result.

Images