Spatial continuous visibility prediction method and device based on multi-modal data analysis

Abstract

The invention provides a spatial continuous visibility prediction method and device based on multi-modal data analysis, and the method comprises the steps: obtaining the multi-modal observation data of a target region, carrying out the image feature coding of the image observation data, generating an atmospheric optical thickness parameter, carrying out the path feature coding of the laser radar scanning data, and carrying out the path feature coding of the laser radar scanning data; generating an extinction coefficient distribution parameter; inputting the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into a multi-modal fusion model to perform optical thickness correction on the atmospheric optical thickness parameter in the horizontal direction to generate a horizontal optical thickness correction parameter, and performing spatial continuous fusion on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter to obtain a multi-modal fusion model; and generating a spatial continuous visibility distribution characteristic of the target area, then converting the spatial continuous visibility distribution characteristic into a visibility prediction result under a polar coordinate system, and outputting a coordinate distribution diagram of the visibility prediction result. According to the invention, the visibility prediction accuracy can be improved.

Description

Technical Field

[0001] The present invention relates to the field of data processing, in particular 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. To obtain spatially continuous visibility, generally, the data of a group of horizontal scans are directly interpolated to obtain the spatial visibility. However, this method has a large error in direct interpolation at a location 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, a method for predicting spatially continuous visibility based on multi-modal data analysis is provided. 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 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 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.

[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 an atmospheric optical thickness parameter, and performing path feature encoding on the lidar scan data to generate an extinction coefficient distribution parameter; a parameter correction module for inputting the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into a multi-modal fusion model, and performing horizontal optical thickness correction on the atmospheric optical thickness parameter through the multi-modal fusion model to generate a horizontal optical thickness correction parameter; a parameter fusion module for performing spatial continuity fusion based on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter to generate a spatial continuous visibility distribution feature of the target area; a visibility prediction module for converting the spatial 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.

[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 extracts the atmospheric optical thickness parameter and the aerosol extinction coefficient distribution parameter by fusing image observation data and lidar scan data respectively, and performs horizontal optical thickness correction and spatial continuity fusion on the two types of parameters based on a multi-modal fusion model, and finally generates a 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 horizontal angular resolution and the lidar in 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 spatially continuously 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 1FIG. shows a schematic diagram of an application scenario provided according to an embodiment of the present invention. 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 includes at least an image acquisition device (such as a pan-tilt-zoom 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 further include one or more databases 130. In certain 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 databases 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 multi-modal data analysis provided by the embodiment of the present invention may specifically include the following steps: Step S100: Obtain 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.

[0013] 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 urban road network, the area surrounding the airport runway, or the continuous space where meteorological monitoring stations are located. Multimodal observation data is a collection of heterogeneous data acquired through different sensors or observation means. In this step, it specifically includes two types: image observation data and lidar scan data. The image observation data is captured by visible light imaging devices (such as pan-tilt cameras), containing two-dimensional image information of the target area in the visible light band, along with additional shooting time and spatial position information. The shooting time is used to identify the exact 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 the 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 height. 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 multimodal observation data set of the target area can be constructed, providing spatially aligned inputs for subsequent feature extraction and fusion.

[0014] 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.

[0015] 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 of visible light by the atmosphere. Its value is determined by the concentration and distribution of components such as aerosols and water vapor, and is specifically 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 of 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 of 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 altitude angle on the optical thickness calculation, the solar altitude 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 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 value is mapped to the extinction coefficient distribution parameter matching the fan-shaped observation sub-region 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.

[0016] Step S300: Input the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into the multi-modal fusion model. The multi-modal fusion model corrects the optical thickness in the horizontal direction for the atmospheric optical thickness parameter to generate a horizontal optical thickness correction parameter.

[0017] The multi-modal fusion model is, for example, a feature alignment and cross-modal association model designed based on a deep learning architecture, which is used to integrate feature parameters from different modalities and eliminate the inconsistencies between modalities. The optical thickness correction in the horizontal direction 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 in the vertical direction of the atmosphere, 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). The correlation weight between the first feature vector and the second feature vector is calculated through a 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 correlation weight to generate a fused feature vector. After being mapped by a fully connected layer, the fused feature vector outputs a horizontal optical thickness correction parameter, which, while retaining the vertical extinction information of the atmosphere, integrates the high-resolution extinction coefficient data of the lidar horizontal scan path, thereby significantly improving the optical thickness estimation accuracy 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.

[0018] 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.

[0019] Spatial continuity fusion is to transform discrete characteristic 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 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 represent the continuous trend of visibility changing with azimuth angle and distance.

[0020] 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.

[0021] 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 element stores the visibility value corresponding to the azimuth angle and distance. The coordinate distribution map visualizes this matrix as a polar coordinate heat map, which is convenient for intuitively displaying 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 the 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 with no 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 for weighted averaging 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.

[0022] As an implementation manner, in step S100, obtaining multi-modal observation data of the target area may specifically include: 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 the absolute value of the timestamp difference less than the preset synchronization threshold as an associated data set with time alignment.

[0023] 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 certain 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, the time-synchronized multi-modal data is screened to ensure the physical consistency of subsequent feature fusion and avoid spatial feature mismatch caused by time misalignment.

[0024] 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.

[0025] 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 circumferential 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, and each cell is uniquely determined by the fan-shaped angle range and the radial distance interval, providing a structured framework for subsequent multi-modal data spatial binding.

[0026] 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 with the corresponding radial distance interval.

[0027] The spatial position 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 central pixel of the image corresponds to the due east direction, the left pixel corresponds to the northward direction, and the right pixel corresponds 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 area 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 azimuth angle of the imaging device (90 degrees) plus the fan angle (30 - 45 degrees), that is, it actually covers 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 area that overlaps with this fan-shaped sub-region (such as a specific rectangular area 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 that the image data is strictly aligned with the lidar scan data in the spatial dimension. For example, if a local observation image block is bound to the "45 - 60 degree, 1000 - 2000 meter" unit, then all its pixels represent the visible light information within this fan angle and distance interval.

[0028] 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 to the local observation image block.

[0029] 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 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). 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), all scan data within the sector angle range are traversed, and the extinction coefficient sampling points with distance values between 500 meters and 1500 meters are selected to form a subset of extinction coefficients. The data in this subset are arranged 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.

[0030] 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.

[0031] 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, its bound local observation image block and subset of extinction coefficients are extracted and arranged 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 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 is strictly spatially aligned 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.

[0032] As an implementation manner, in step S200, image feature encoding is performed on the image observation data to generate an atmospheric optical thickness parameter, which may specifically include: Step S210: Extract the pixel brightness distribution sequence of the sky region in the visible light image, where the pixel brightness distribution sequence includes the brightness means of the sky region at different azimuth angles.

[0033] 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, and is usually a set of pixels at the top of the image or 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 pixels of the sky region 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.

[0034] 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.

[0035] 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 (such as 0 - 15 degrees, 15 - 30 degrees, etc.) in the pixel brightness distribution sequence, 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 as the azimuth angle increases, that is, the atmospheric extinction effect increases; a positive gradient value indicates that the brightness increases as the azimuth angle increases, that is, the atmospheric transparency improves. By traversing all adjacent azimuth angle partitions, an optical attenuation gradient distribution sequence covering the complete 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 non-uniformity. For example, the presence of consecutive negative gradients within a certain azimuth angle range indicates the existence of an aerosol aggregation area in that direction.

[0036] 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.

[0037] 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 a high aerosol concentration in that direction and output a relatively high optical thickness value (such as 0.8). During the training phase, 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 and 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.

[0038] 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.

[0039] 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.

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

[0041] Spatial projection is to map the azimuth angle dimension parameters of the corrected optical thickness distribution to the sector 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 sector 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 sector 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 sector 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 sector 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 sector sub-region" corresponds to the optical thickness parameter 0.566.

[0042] 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: 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.

[0043] 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 per 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 segmented 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, its 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.

[0044] 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.

[0045] 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 the 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.

[0046] 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.

[0047] 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.

[0048] 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.

[0049] 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 changing with distance. For example, the contribution value in the short-distance segment is relatively low (low aerosol concentration), and the contribution value gradually increases in the long-distance segment (aerosol concentration accumulation).

[0050] 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 matching the fan-shaped observation sub-region.

[0051] 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 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 into 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, 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 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.

[0052] As an implementation manner, in step S300, the multi-modal fusion model is used to correct the optical thickness in the horizontal direction of the atmospheric optical thickness parameter to generate a horizontal optical thickness correction parameter, which may specifically include: 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.

[0053] 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 in a range of azimuth angles (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 in 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 uses 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 uses 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.

[0054] Step S320: Calculate the association 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 association weight to generate a fused feature vector.

[0055] The cross-attention mechanism is the core module of the multi-modal fusion model. Its function is to model the spatial correlation between different modal features 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 association weight between the query vector and the key vector through matrix operations. For example, for a certain fan-shaped sub-region feature in the first feature vector (such as the atmospheric optical thickness parameter in the 60-75 degree fan-shaped 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 association weight through the Softmax function. The association weight reflects the spatial dependence strength between the atmospheric optical thickness parameter of this fan-shaped 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 fan-shaped sub-region has a high association 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 fan-shaped sub-region and the horizontal extinction characteristics in the 2000-2050 meter distance segment. Based on the association weight, a weighted sum of the second feature vector is performed to generate a context-enhanced feature for the first feature vector, and then it is connected with the original first feature vector through a residual connection 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.

[0056] 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.

[0057] 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 sector 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 a 24-dimensional output corresponding to the correction parameters of 24 sector sub-regions), and then performs 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 sector 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, which 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 sector 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 sector 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 accurate 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 sector 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.

[0058] As an implementation manner, the training process of the above multi-modal fusion model may include the following steps: 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.

[0059] 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 set of visible light images of the target area collected by a visible light imaging device within a certain past time period. Each image is attached with an accurate capture timestamp (such as "2022-05-10 14:30:00 UTC") and spatial location information (such as latitude and longitude coordinates, capture azimuth angle, and pitch angle parameters). The historical lidar scan data is the aerosol extinction coefficient distribution data collected by a lidar device along a horizontal scan path during the same period. Each scan data records the scan start timestamp, the horizontal scan angle range (such as 0-360 degrees), and the sequence of extinction coefficient sampling points along the radial path. The visibility ground truth distribution map is generated by professional meteorological observation equipment (such as a transmissometer or a multi-point scattering sensor network) under the same spatio-temporal conditions. It 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 a Cartesian coordinate system or a polar coordinate system. For example, the historical image observational data may include a visible light image captured at 14:30 on May 10, 2022. The corresponding lidar scan data is the extinction coefficient distribution of a 0-360 degree horizontal scan at the same moment, and the visibility ground truth 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 visibility ground truth distribution map) to support the accuracy of the supervision signal for model training.

[0060] 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.

[0061] The feature extraction process of historical image observation data is the same as steps S200 - S250. 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 means 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 for 24 sector - shaped sub - regions is generated, and each parameter represents the vertical extinction integral value in the corresponding azimuth range (such as 30 - 45 degrees). The feature extraction process of historical lidar scan data is the same as 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 for 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.

[0062] 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 the predicted visibility distribution features.

[0063] 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 to generate 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, when 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 500 - 1000 - meter distance segment in this region.

[0064] 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.

[0065] 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, then 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 grid cells to obtain the mean squared error (MSE) as the spatial consistency loss value. The model training uses the backpropagation algorithm to calculate the gradient based on the loss value and update the network weights (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 consecutive multiple training epochs, which is regarded as convergence.

[0066] As an implementation manner, in Step S40, calculating the spatial consistency loss value between the predicted visibility distribution feature and the true visibility distribution map may specifically include: 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.

[0067] 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.

[0068] 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: 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 angle resolution. The polar coordinate grid parameters include the angle dimension division sequence and the distance dimension division sequence.

[0069] The spatial coverage range is determined by the maximum radial distance (such as 5000 m) and the angle coverage range (such as 0 - 360 degrees) of the predicted visibility distribution characteristics. The preset angle resolution determines the division interval of the angle 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 angle resolution of 15 degrees generates 24 angle partitions. The polar coordinate grid parameters are defined as the angle dimension division sequence [0°, 15°, 30°…345°] and the distance dimension division sequence [0 m, 50 m, 100 m…5000 m].

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

[0071] The conversion formula from Cartesian coordinates to polar coordinates is: angle 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 angle 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 45 - 60 - degree angular partition and the 2950 - 3000 - m distance segment.

[0072] 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.

[0073] A blank cell refers to a cell in the polar coordinate grid without a predicted visibility value mapped, which may be caused by the sparsity of the original prediction data or coordinate conversion errors. The 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 in 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.

[0074] Step S414: Synchronously convert each true visibility value in the true visibility distribution map to 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.

[0075] 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 45 - 60 - degree angular partition and the 2100 - 2150 - m distance segment. If the true visibility of this cell is 8500 m, the corresponding grid cell of the second visibility matrix is assigned 8500 m.

[0076] 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.

[0077] 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 in the 30 - 45 degree angular partition of the true value distribution, 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 the same spatial structure as the first visibility matrix, thus supporting accurate loss value calculation.

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

[0079] 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 (such as 24 angular partitions × 100 distance segments = 2400 cells), and accumulate all differences to generate the total loss value.

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

[0081] 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 areas. 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.

[0082] 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 characteristics of the target area, which may specifically include: 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.

[0083] 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 scanning 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). As the initial input, the original resolution of the horizontal optical thickness correction parameter may not be consistent with the preset angular resolution (for example, 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 are divided with high resolution in the angular and distance dimensions, providing refined input for subsequent spatial continuity fusion.

[0084] 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.

[0085] The angle - distance grid is formed by the intersection of the angle dimension (0 - 360 degrees, at 1 - degree intervals) and the distance dimension (0 - 5000 meters, at 10 - meter intervals). 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 angles 5 - 6 degrees and distance 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, 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 at 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 need to be further optimized.

[0086] 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.

[0087] 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 angle 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 the context information in a larger range, and finally the output layer maps the feature map to the visibility value space through a 1×1 convolution. During the training process, the spatial continuity optimization network learns to suppress the high-frequency noise in the interpolation results (such as the abnormally high value of a single grid) and enhance the spatial continuity (such as the gradual change of visibility between adjacent angles and distances). For example, if the visibility of the cell at an angle of 5 degrees and a distance of 100 meters in the interpolation result is 0.061, while the visibility of the adjacent cell at an angle of 5 degrees and a distance of 110 meters is 0.12 (a sudden change caused by 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.

[0088] Among them, the training process of the above spatial continuity optimization network may include the following steps: Step S4301: Obtain synthetic visibility training data, which includes simulated visibility interpolation results and corresponding continuous visibility distribution labels.

[0089] The synthetic visibility training data is generated through numerical simulation or adding noise to real data. The numerical simulation method generates idealized continuous visibility distribution labels based on an atmospheric physical 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, in a continuous visibility distribution label, the visibility in a certain area decreases linearly from 5000 meters (angle 0 - 90 degrees) to 1000 meters (angle 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 a variety of atmospheric conditions (such as fog, haze, clear sky) and noise patterns to ensure the generalization ability of the model.

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

[0091] The initial spatial continuity optimization network is an untrained CNN model with the same structure as 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 the convolutional kernel. 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.

[0092] 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.

[0093] 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 the 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.

[0094] As an implementation, in step S500, converting the spatially continuous visibility distribution feature into a visibility prediction result in the polar coordinate system may specifically include: 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.

[0095] 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 circumference 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 with a 50-meter interval. 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.

[0096] 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.

[0097] 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 transformation 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 azimuth. Through this step, the Cartesian coordinates of all 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.

[0098] 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.

[0099] 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 in step S520 into the corresponding grid cells respectively: the angular component θ = 53.13 degrees belongs to the 45 - 60 degree angular axis grid cell, and the distance component r = 5000 meters belongs to the 4950 - 5000 meter distance axis grid cell. 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 45 - 60 degree, 4950 - 5000 meter cell), the mean 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 regions in the initial mapping table, which need to be filled by interpolation.

[0100] 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 adjacent interpolation to generate the interpolated visibility value of the current grid cell.

[0101] Two - way adjacent 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 30 - 45 degrees in angle and 2000 - 2050 meters in distance), extract the visibility values of adjacent angular cells (such as 15 - 30 degrees and 45 - 60 degrees) along the angular axis direction, and extract the visibility values of 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 adjacent angular cells are 8000 meters (15 - 30 degrees) and 8200 meters (45 - 60 degrees), and the visibility values of 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 of the angular - direction mean value (8000 + 8200) / 2 = 8100 meters and the distance - direction mean value (7900 + 8100) / 2 = 8000 meters, that is, 8050 meters. If there is no valid data in one 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.

[0102] 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 entire polar coordinate range.

[0103] 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 according to the sequences divided by the angle axis (0°, 15°…345°) and the distance axis (0m, 50m…5000m). For example, the first column of grid cells corresponding to the 0° on the angle axis stores the visibility values of the distance segments from 0 - 50 meters, 50 - 100 meters…4950 - 5000 meters in sequence. The second column corresponding to the 15° on the angle axis stores the visibility values of the same distance segments, until the last column corresponding to the 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 the contour interval represents the level of the visibility value (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 of the continuous distribution in the entire space of the target area. For example, it can identify a fog area with a 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.

[0104] In a possible implementation, after outputting the coordinate distribution map of the visibility prediction results, it may further include the process of generating dynamic visibility warnings, which may specifically include the following steps: 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.

[0105] 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 the historical visibility prediction data. For example, it is collected once every 5 minutes, and 12 groups 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.

[0106] 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.

[0107] 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 angle axis units × number of distance axis units.

[0108] 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.

[0109] 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, if the fluctuation amplitude of the cell at an angle of 30 - 45 degrees and a distance of 2000 - 2050 m at time point T2 is -600 m (lower than -500 m), 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.

[0110] 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.

[0111] Spatial continuity clustering uses a regional growing algorithm to merge adjacent abnormal fluctuation units to form a continuous region. The definition of proximity includes proximity in angle and distance dimensions: at the same time point, if the angle difference between two abnormal units is less than the preset angle tolerance (such as 15 degrees) and the distance difference is less than the 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 the angle difference is 15 degrees and the distance difference is 50 meters, which meets the proximity conditions and is clustered into the same abnormal fluctuation area. The regional boundary is generated by extracting the outer contour of the clustering area, for example, using the convex hull algorithm to determine the minimum enclosing polygon. Based on the original coordinate distribution map, the dynamic warning annotation map marks the spatial position and time range of the abnormal fluctuation area with highlighted color blocks or contour lines. For example, a red translucent overlay is used to indicate that there is abnormally low visibility fluctuation in the area with an angle of 30-90 degrees and a distance of 1500-3000 meters during the period from 14:05 to 14:20.

[0112] In a possible implementation manner, after the dynamic warning annotation map is generated, the embodiment of the present invention may further include a process of verifying and correcting the abnormal area, which may specifically include the following steps: Step S1000: Obtain the historical visibility data corresponding to the abnormal fluctuation area, and extract the visibility mean sequence of the same spatial position in the historical period.

[0113] The historical visibility data for the same period refers to the meteorological observation data for the same period (such as the same month or time period) in the past few years. For example, if the current date is 14:00-15:00 on November 15, 2023, the historical data for the same period are the measured visibility values ​​from 14:00 to 15:00 on November 15 every year from 2020 to 2022. For each grid cell in the abnormal fluctuation area (such as angle 30-45 degrees, distance 2000-2050 meters), the visibility values ​​of the same spatial position in the same period in history are extracted, and their mean sequence is calculated. For example, the visibility of the cell is 8500 meters in 2020, 8200 meters in 2021, and 8000 meters in 2022. The mean value for the same period in history is (8500+8200+8000) / 3≈8233 meters. The mean sequence reflects the typical visibility level of the location in the same period in history, providing a benchmark for deviation analysis.

[0114] Step S1100: compare the real-time visibility value with the visibility mean sequence grid unit by grid unit to generate a visibility deviation distribution map.

[0115] The formula for calculating the deviation degree is: Deviation degree = (Real-time visibility value - Mean value of the same period in history) / Mean value of the same period in history × 100%. For example, if the real-time visibility is 7500 meters and the mean value of the same period in history 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 level of the same period in history, 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%).

[0116] 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.

[0117] 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 mean value of the same period in history), 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% three times in a row. This dual-threshold mechanism can filter out short-term fluctuations and focus on long-term anomalies.

[0118] 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.

[0119] 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, distinguishes 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.

[0120] In a possible implementation manner, after generating the graded early warning distribution map, the embodiments of the present invention may further include the process of real-time visibility update and feedback, which may specifically include the following steps: 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.

[0121] 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, the visibility in the high-confidence anomaly area at time point T1 (14:00) is 8000 meters, at T2 (14:05) it is 7500 meters, and the time interval is 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.

[0122] 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.

[0123] 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.

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

[0125] 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 on the graded early warning distribution map.

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

[0127] 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 (the intersection of 30 - 60 degrees and 30 - 90 degrees) × (the intersection of 2000 - 2500 m and 1500 - 3000 m) / the area of the 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.

[0128] 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 constituent structures: 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; A feature encoding module 920, configured to perform image feature encoding on the image observation data to generate an atmospheric optical thickness parameter, and perform path feature encoding on the lidar scan data to generate an extinction coefficient distribution parameter; A parameter correction module 930, configured to input the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into a multi-modal fusion model, and perform horizontal optical thickness correction on the atmospheric optical thickness parameter through the multi-modal fusion model to generate a horizontal optical thickness correction parameter; A parameter fusion module 940, configured to perform spatial continuity fusion based on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter to generate a spatial continuous visibility distribution feature of the target area; A visibility prediction module 950, configured to convert the spatial continuous visibility distribution feature into a visibility prediction result in the polar coordinate system and output a coordinate distribution map of the visibility prediction result.

[0129] 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 spatial continuous visibility prediction method based on multimodal data analysis, characterized in that: The method comprises: Acquire multimodal observation data of the target area, the multimodal observation data including image observation data and laser radar scanning data, wherein the image observation data includes a visible light image of the target area and corresponding shooting time and spatial position information, and the laser radar scanning data includes aerosol extinction coefficient distribution data on a horizontal scanning 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 laser radar scanning data to generate extinction coefficient distribution parameters; Inputting the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into a multimodal fusion model, performing optical thickness correction in a horizontal direction on the atmospheric optical thickness parameter through the multimodal fusion model, and generating a horizontal optical thickness correction parameter; Based on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter, spatial continuity fusion is performed to generate spatial continuous visibility distribution characteristics of the target area; The spatial continuous visibility distribution characteristics are converted into visibility prediction results in a polar coordinate system, and a coordinate distribution diagram of the visibility prediction results is output.

2. The method according to claim 1, characterized in that The obtaining of multimodal observation data of the target area includes: Extracting the shooting timestamp of the image observation data and the scanning timestamp of the laser radar scanning data, and marking the image observation data and the laser radar scanning data whose timestamp difference is less than a preset synchronization threshold as a time-aligned associated data set; Based on the horizontal scanning angle range of the laser radar scanning data in the associated data set, the target area is divided into a plurality of fan-shaped observation sub-areas at preset angle intervals, and a corresponding radial distance interval is allocated to each of the fan-shaped observation sub-areas; According to the spatial position information of the image observation data in the associated data set, a pixel area in the visible light image that overlaps with the angle range of the fan-shaped observation sub-area is extracted as a local observation image block, and the local observation image block is spatially bound to the corresponding radial distance interval; Performing range interval filtering on the laser radar scanning data in each of the fan-shaped observation sub-areas, retaining the aerosol extinction coefficient distribution data in the radial distance interval, and generating an extinction coefficient subset corresponding to the local observation image block one by one; The local observation image blocks and the extinction coefficient subsets are arranged in the angular order of the fan-shaped observation sub-regions to generate a spatially aligned multimodal observation data sequence.

3. The method according to claim 1, characterized in that The step of encoding the image features of the image observation data to generate atmospheric optical thickness parameters includes: Extracting a pixel brightness distribution sequence of a sky region in the visible light image, wherein the pixel brightness distribution sequence includes brightness means of the sky region at different azimuth angles; Generating an optical attenuation gradient distribution of the visible light image in a horizontal direction according to a variation trend of brightness mean values ​​of adjacent azimuth angles in the pixel brightness distribution sequence; Inputting the optical attenuation gradient distribution into a pre-trained atmospheric scattering correction network, mapping the nonlinear relationship between the optical attenuation gradient distribution and the atmospheric scattering parameters through a multi-layer perception structure in the atmospheric scattering correction network, and generating an initial atmospheric optical thickness distribution; Acquire the solar altitude angle parameter corresponding to the shooting time, perform angle compensation on the initial atmospheric optical thickness distribution through the time correction module in the atmospheric scattering correction network, and generate a corrected optical thickness distribution aligned with the solar position; Performing spatial projection on the corrected optical thickness distribution on the fan-shaped observation sub-area of ​​the target area to generate local atmospheric optical thickness parameters corresponding to each fan-shaped observation sub-area; The step of encoding the path characteristics of the laser radar scanning data to generate extinction coefficient distribution parameters includes: Performing scanning path segmentation on the aerosol extinction coefficient distribution data on the horizontal scanning path according to a preset distance resolution to generate a plurality of radial distance segments corresponding to the fan-shaped observation sub-areas; For each radial distance segment, extracting continuous extinction coefficient sampling points located in the distance segment from the aerosol extinction coefficient distribution data, and calculating the mean extinction coefficient of the continuous extinction coefficient sampling points; Generate an extinction integral contribution value of each radial distance segment according to the mean extinction coefficient and the length of the corresponding radial distance segment; Arranging the extinction integral contribution values ​​according to the spatial order of the radial distance segments to generate an extinction integral value sequence on the horizontal scanning path; Based on the conversion relationship between the extinction integral value sequence and the preset visibility, the extinction integral value sequence is mapped into an extinction coefficient distribution parameter matching the fan-shaped observation sub-area.

4. The method according to claim 1, characterized in that: The method of performing horizontal optical thickness correction on the atmospheric optical thickness parameter by using the multimodal fusion model to generate a horizontal optical thickness correction parameter includes: Inputting the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into the feature alignment layer of the multimodal fusion model to generate a first eigenvector of the atmospheric optical thickness parameter in the horizontal angle and a second eigenvector of the extinction coefficient distribution parameter in the distance dimension; By using the cross attention mechanism of the multimodal fusion model, the association weight between the first feature vector and the second feature vector is calculated, and the first feature vector is weightedly fused based on the association weight to generate a fused feature vector; The fused feature vector is input into the fully connected layer of the multimodal 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 multimodal fusion model includes: Acquire a historical multimodal training data set, wherein the historical multimodal training data set includes historical image observation data, historical lidar scanning data, and a corresponding true visibility distribution map; Performing feature extraction on the historical image observation data to generate historical atmospheric optical thickness parameters, and performing feature extraction on the historical lidar scanning data to generate historical extinction coefficient distribution parameters; Inputting the historical atmospheric optical thickness parameter and the historical extinction coefficient distribution parameter into the initial multimodal fusion model to generate predicted visibility distribution characteristics; A spatial consistency loss value between the predicted visibility distribution feature and the true visibility distribution map is calculated, and parameters of the initial multimodal fusion model are adjusted 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 calculating of the spatial consistency loss value between the predicted visibility distribution feature and the true visibility distribution map comprises: Converting the predicted visibility distribution characteristics and the true visibility distribution map into a first visibility matrix and a second visibility matrix in a polar coordinate system respectively; Calculating the visibility difference between the first visibility matrix and the second visibility matrix in the same angle and distance unit; The visibility difference is normalized to generate the spatial consistency loss value.

7. The method according to claim 6, characterized in that The step of converting the predicted visibility distribution feature and the visibility true value distribution map into a first visibility matrix and a second visibility matrix in a polar coordinate system respectively includes: Extracting the spatial coverage of the predicted visibility distribution feature, and generating polar coordinate grid parameters based on the maximum radial distance of the spatial coverage and a preset angular resolution, wherein 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 feature into an angle value and a distance value in a polar coordinate system according to its spatial coordinates in a Cartesian coordinate system, and map the angle value and the distance value to corresponding polar coordinate grid cells according to the angle dimension partition sequence and the distance dimension partition sequence; For blank cells in the polar coordinate grid cells that are not covered by the predicted visibility distribution characteristics, radial interpolation filling is performed based on visibility values ​​of adjacent grid cells to generate a first visibility matrix covering the complete polar coordinate grid; Synchronously converting each true visibility value in the true visibility distribution map into a true angle value and a true distance value in a polar coordinate system according to its spatial coordinates, and mapping the true angle value and the true distance value to corresponding polar coordinate grid cells based on the same polar coordinate grid parameters; The blank cells in the polar coordinate grid that are not covered by the true visibility values ​​are filled 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, characterized in that The spatially continuous fusion based on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter to generate the spatially continuous visibility distribution characteristics of the target area includes: Dividing the horizontal optical thickness correction parameter into a plurality of angle sub-parameters according to a preset angle resolution, and dividing the extinction coefficient distribution parameter into a plurality of distance sub-parameters according to a preset distance resolution; Performing linear interpolation on each of the angle sub-parameters and the corresponding distance sub-parameters to generate visibility interpolation results within a continuous angle-distance grid; Inputting the visibility interpolation result into a spatial continuity optimization network, smoothing the visibility interpolation results of adjacent grids through the spatial continuity optimization network, and generating the spatial continuous visibility distribution feature; Wherein, the training method of the spatial continuity optimization network includes: Acquire synthetic visibility training data, where the synthetic visibility training data includes simulated visibility interpolation results and corresponding continuous visibility distribution labels; Inputting the simulated visibility interpolation result into the initial spatial continuity optimization network to generate a predicted continuous visibility distribution; A gradient difference loss value between the predicted continuous visibility distribution and the continuous visibility distribution label is calculated, and parameters of the initial spatial continuity optimization network are adjusted based on the gradient difference loss value.

9. The method according to claim 2 or 3, characterized in that: The step of converting the spatial continuous visibility distribution characteristics into a visibility prediction result in a polar coordinate system includes: Based on the maximum radial distance of the spatial continuous visibility distribution feature and the preset angle interval of the fan-shaped observation sub-area, generating polar coordinate grid division parameters, the polar coordinate grid division parameters including an angle axis division sequence and a distance axis division sequence; Extracting the original coordinates of each visibility value in the spatial continuous visibility distribution feature in a Cartesian coordinate system, and converting the original coordinates into an angle component and a distance component in a polar coordinate system; Mapping the angle component to the corresponding angle axis grid unit according to the angle axis division sequence, and mapping the distance component to the corresponding distance axis grid unit according to the distance axis division sequence, to generate an initial polar coordinate visibility distribution mapping table; Traversing each grid cell in the initial polar coordinate visibility distribution mapping table, if there is no mapped visibility value in the current grid cell, extracting the visibility values ​​in the adjacent angle direction and distance direction of the current grid cell for bidirectional neighbor interpolation to generate an interpolated visibility value of the current grid cell; The visibility values ​​of all grid cells after interpolation are rearranged according to the angle axis division sequence and the distance axis division sequence to generate a coordinate distribution diagram of the visibility prediction results covering the complete polar coordinate range.

10. A spatial continuous visibility prediction device based on multimodal data analysis, characterized in that: include: A data acquisition module, used to acquire multimodal observation data of a target area, wherein the multimodal observation data includes image observation data and laser radar scanning data, wherein the image observation data includes a visible light image of the target area and corresponding shooting time and spatial position information, and the laser radar scanning data includes aerosol extinction coefficient distribution data on a horizontal scanning path of the target area; A feature encoding module, used to perform image feature encoding on the image observation data to generate atmospheric optical thickness parameters, and to perform path feature encoding on the laser radar scanning data to generate extinction coefficient distribution parameters; A parameter correction module, used for inputting the atmospheric optical thickness parameter and the extinction coefficient distribution parameter into a multimodal fusion model, performing optical thickness correction in a horizontal direction on the atmospheric optical thickness parameter through the multimodal fusion model, and generating a horizontal optical thickness correction parameter; A parameter fusion module, used for performing spatial continuity fusion based on the horizontal optical thickness correction parameter and the extinction coefficient distribution parameter to generate a spatial continuous visibility distribution feature of the target area; The visibility prediction module is used to convert the spatial continuous visibility distribution characteristics into visibility prediction results in a polar coordinate system and output a coordinate distribution diagram of the visibility prediction results.

Images