A tunnel entrance black hole effect road condition image monitoring method

Abstract

The application provides a tunnel entrance black hole effect road condition image monitoring method, relates to the technical field of visual intelligent detection and image processing, and acquires a three-dimensional prior tensor representing the physical law of macroscopic radiation attenuation; the original image is desensitized and a dynamic photometric trajectory differential quotient set based on speed normalization is extracted; the prior tensor is discretized to the space-time topology space of the differential quotient set, a residual error matrix is generated by performing a Bayesian residual assimilation calculation; a forward detection signal is generated based on the dynamic allocation of a nonlinear fusion weight based on the residual error matrix, and a high-frequency variation component is extracted as a regularization penalty term to trigger reverse iteration of a local climate parameter of the model. The application overcomes the limitation that traditional static two-dimensional visual observation is easily blinded by environmental interference, establishes a macro-micro coordinated two-way evolution closed loop, and realizes high-fidelity brightness detection under meteorological and dynamic interference.

Description

Technical Field

[0001] This invention relates to the field of visual intelligent detection and image processing technology, specifically a method for monitoring road conditions using images showing the black hole effect at tunnel entrances. Background Technology

[0002] With the deepening of the digital transformation of transportation infrastructure, intelligent perception systems based on computer vision and pattern recognition technologies are playing an increasingly crucial role in complex road network environments. Accurately acquiring and analyzing road condition image data in confined spaces with highly dynamic lighting characteristics, such as tunnel entrances, has irreplaceable engineering value for preventing the "black hole effect" and ensuring driving safety. This type of visual intelligent detection not only requires processing massive amounts of high-dimensional image data but also demands that the perception algorithm possess high stability and physical self-consistency under non-ideal conditions.

[0003] Current mainstream methods for detecting the "black hole effect" at tunnel entrances primarily rely on static analysis of relative grayscale in two-dimensional plane images, which limits their adaptability to complex scenarios. For example, Chinese patent application CN109300160A discloses an intelligent, rapid method for detecting the tunnel black hole effect by acquiring conversion coefficients and calculating absolute brightness to determine the effect. While this technology improves detection efficiency under ideal lighting conditions, it still faces two core technical challenges: at the microscopic observation level, it heavily relies on offline preset static brightness conversion coefficients, making it prone to baseline drift and data contamination due to interference in complex dynamic road conditions (vehicle obstruction, strong backlighting and glare); at the macroscopic prediction level, it lacks constraints from three-dimensional spatial physical laws, making it highly susceptible to complete blindness when extreme weather (heavy rain, localized dark fog) causes severe attenuation or loss of visual features, rendering a single visual perception algorithm prone to complete blindness. These limitations collectively point to a clear technical gap: existing solutions struggle to establish a steady-state physical connection between microscopic dynamic visual interference and abrupt changes in the macroscopic physical environment. Summary of the Invention

[0004] The purpose of this invention is to provide a method for monitoring road condition images of the "black hole effect" at tunnel entrances, which decouples and re-integrates macroscopic meteorological knowledge graphs with microscopic dynamic photometric trajectories in time and space. By introducing a three-dimensional prior radiation tensor constrained by Maxwell's equations and combining it with Bayesian residual assimilation calculations, this invention establishes a bidirectional collaborative evolution mechanism of "observation optimization based on prior prediction and self-evolution of prior models based on microscopic observations," thereby reconstructing the mapping relationship between visual features and physical laws at the underlying logic level to solve the problems mentioned in the background art.

[0005] To achieve the above objectives, the present invention provides the following technical solution:

[0006] A method for monitoring road condition images of the black hole effect at tunnel entrances, comprising the following steps:

[0007] S1: Obtain the three-dimensional radiation attenuation prior tensor that represents the macroscopic radiation attenuation physical law within the tunnel entrance transition space defined by the target geographic coordinates and local terrain structure.

[0008] S2: Obtain the original image sequence characterizing the spatiotemporal changes of photometric surface of dynamic target within the specific spatial environment; perform irreversible edge anonymization and desensitization on the original image sequence to remove identity features; and extract the set of dynamic photometric trajectory differential quotients based on the desensitized image sequence.

[0009] S3: Discretize the three-dimensional radiation attenuation prior tensor into the spatiotemporal topological space corresponding to the set of differential quotients of dynamic photometric trajectories, and perform Bayesian residual assimilation calculation to generate a spatiotemporal residual distribution matrix that maps the difference between the three-dimensional radiation attenuation prior tensor and the set of differential quotients of dynamic photometric trajectories.

[0010] S4: Based on the spatiotemporal residual distribution matrix, allocate the fusion weight parameters of the three-dimensional radiation attenuation prior tensor and the set of differential quotients of the dynamic photometric trajectory, and perform a weighted fusion operation based on the fusion weight parameters to generate a positive detection indication signal characterizing the global adaptive calibration absolute brightness attenuation distribution.

[0011] S5: Extract the high-frequency variation components in the spatiotemporal residual distribution matrix, and use the high-frequency variation components as regularization penalty terms to generate a reverse evolution indication signal for triggering the iterative update of the local climate parameters of the three-dimensional radiation attenuation prior tensor in the next calculation cycle.

[0012] Compared with the prior art, the beneficial effects of the present invention are:

[0013] This invention performs irreversible edge anonymization desensitization on the original image sequence and extracts the set of differential quotients of dynamic photometric trajectories (S2). It then introduces a high light repulsion threshold and a transient velocity scalar for dimensional normalization. This technique effectively filters out invalid specular reflections and eliminates the interference of target motion speed differences on photometric attenuation calibration. It fundamentally overcomes the baseline drift problem caused by strong backlighting and complex traffic flow, which is common in traditional methods. This improves the physical consistency of the microscopic probe, thereby enhancing the robustness and data purity of microscopic dynamic observations.

[0014] This invention constructs a macroscopic logical knowledge graph containing ephemeris and meteorological features, and performs graph neural network topological reasoning (S1) based on simplified boundary condition constraints according to Maxwell's radiative transfer equation to calculate a three-dimensional radiative attenuation prior tensor with spatial depth attributes. This method provides three-dimensional physical constraints for visual detection, ensuring an absolute theoretical calculation benchmark even under extreme conditions where basic visual features are lost due to heavy rain or dark fog, maintaining the integrity of environmental perception information entropy, and overcoming the visual blindness defects under extreme weather conditions.

[0015] This invention generates a spatiotemporal residual distribution matrix by performing Bayesian residual assimilation calculation (S3), and allocates fusion weights based on a nonlinear relational mapping function constructed according to the total number of effective samples and the global residual norm (S4). When the number of probes is lower than the detection threshold, this mechanism can dynamically adjust the weight transition based on the steepness of environmental disturbances, avoiding the step oscillations caused by traditional Boolean hard switching logic, and providing a highly sensitive and absolutely smooth safety control degradation curve to achieve cross-scale smooth degradation and high stability of system control;

[0016] This invention extracts high-frequency variation components from the spatiotemporal residual distribution matrix and uses them as regularization penalty terms (S5) when the spatial clustering index convergence condition is met. By forcibly updating the local aerosol scattering coefficients in the graph neural network through error backpropagation, this inverse evolution method achieves a leap from pure data-driven to physics-inspired approaches. It can complete the adaptive iteration of model parameters without human intervention, establish robust synergistic gains across evolutionary cycles, and form a synergistic gain closed loop between physical parameters and model parameters. Attached Figure Description

[0017] Figure 1 A hierarchical architecture and core data flow technology roadmap for a method of monitoring road conditions with the black hole effect at tunnel entrances;

[0018] Figure 2 This is a technical roadmap for generating the spatiotemporal residual distribution matrix of this invention;

[0019] Figure 3 This is a technical roadmap for generating the reverse evolution indicator signal of the present invention;

[0020] Figure 4 The diagram shows a numerical comparison and verification of the traditional Boolean scheme and the nonlinear degradation weight allocation logic of this invention. Detailed Implementation

[0021] To make the objectives, technical solutions, and advantages of this invention clearer, the invention will be further described in detail below with reference to the accompanying drawings and embodiments. It should be understood that the specific embodiments described herein are merely illustrative and not intended to limit the invention.

[0022] It is understood that the terms “first,” “second,” etc., used in this application may be used herein to describe various elements, but unless otherwise stated, these elements are not limited by these terms. These terms are used only to distinguish one element from another.

[0023] Example 1:

[0024] Please see Figures 1 to 4 The present invention provides a technical solution:

[0025] A method for monitoring road condition images of the black hole effect at tunnel entrances includes the following steps:

[0026] S1: Obtain the three-dimensional radiation attenuation prior tensor that represents the macroscopic radiation attenuation physical law within the tunnel entrance transition space defined by the target geographic coordinates and local terrain structure.

[0027] S2: Obtain the original image sequence characterizing the spatiotemporal changes of photometric surface of dynamic target within the specific spatial environment; perform irreversible edge anonymization and desensitization on the original image sequence to remove identity features; and extract the set of dynamic photometric trajectory differential quotients based on the desensitized image sequence.

[0028] S3: Discretize the three-dimensional radiation attenuation prior tensor into the spatiotemporal topological space corresponding to the set of differential quotients of dynamic photometric trajectories, and perform Bayesian residual assimilation calculation to generate a spatiotemporal residual distribution matrix that maps the difference between the three-dimensional radiation attenuation prior tensor and the set of differential quotients of dynamic photometric trajectories.

[0029] S4: Based on the spatiotemporal residual distribution matrix, allocate the fusion weight parameters of the three-dimensional radiation attenuation prior tensor and the set of differential quotients of the dynamic photometric trajectory, and perform a weighted fusion operation based on the fusion weight parameters to generate a positive detection indication signal characterizing the global adaptive calibration absolute brightness attenuation distribution.

[0030] S5: Extract the high-frequency variation components in the spatiotemporal residual distribution matrix, and use the high-frequency variation components as regularization penalty terms to generate a reverse evolution indication signal for triggering the iterative update of the local climate parameters of the three-dimensional radiation attenuation prior tensor in the next calculation cycle.

[0031] Further defining step S1, the step of obtaining the three-dimensional radiation attenuation prior tensor further includes:

[0032] Obtain macroscopic logical knowledge graph data containing target geographic coordinates, ephemeris status, and real-time meteorological cloud cover status;

[0033] Using the macroscopic logical knowledge graph data as input features, a graph neural network topological reasoning process based on simplified boundary condition constraints according to Maxwell's radiative transfer equation is executed to calculate the three-dimensional radiative attenuation prior tensor with spatial depth attributes.

[0034] Further defining the steps of executing the graph neural network topological reasoning process based on simplified boundary condition constraints according to Maxwell's radiative transfer equation, the steps further include: extracting the global meteorological feature map and the local terrain structure map from the macroscopic logical knowledge graph data; calculating the cross-scale attention weight matrix between the global meteorological feature map and the local terrain structure map; and performing node feature aggregation on the initial three-dimensional radiation field based on the cross-scale attention weight matrix to generate the three-dimensional radiation attenuation prior tensor.

[0035] Further defining step S2, the step of extracting the set of dynamic photometric trajectory differential quotients based on the desensitized image sequence further includes:

[0036] Extract the dynamic target bounding box sequence from the desensitized continuous temporal images;

[0037] For each frame target region in the dynamic target bounding box sequence, calculate its monochrome variation coefficient in multiple color channels in the first color space;

[0038] An arithmetic mean is performed on the monochrome coefficients of variation across multiple color channels to generate the pixel coefficients of variation for the dynamic target.

[0039] The pixel variation coefficient is compared with a preset specular repulsion threshold; when the pixel variation coefficient is lower than the preset specular repulsion threshold, the apparent gray-scale temporal difference of the same dynamic target on the spatial displacement vector near the tunnel entrance is calculated to generate the set of differential quotients of the dynamic photometric trajectory.

[0040] Further defining the step of calculating the apparent gray-scale temporal difference of the same dynamic target on the spatial displacement vector near the tunnel entrance, the step further includes: extracting the motion optical flow vector of the same dynamic target between consecutive frames; calculating the transient velocity scalar based on the motion optical flow vector; introducing a dynamic normalization parameter, and performing a multiplication operation between the transient velocity scalar and the dynamic normalization parameter to generate a dimensionless velocity penalty factor.

[0041] The ratio of the apparent gray-scale temporal difference to the velocity penalty factor is calculated, and this ratio is mapped to a preset normalization interval to generate the set of differential quotients of the dynamic photometric trajectory. This eliminates the interference of target motion velocity differences on photometric attenuation calibration.

[0042] Further defining step S3, the step of performing Bayesian residual assimilation calculation to generate the spatiotemporal residual distribution matrix of the mapping difference further includes:

[0043] A four-dimensional spatiotemporal interpolation grid is constructed based on a preset camera intrinsic parameter transformation matrix; the three-dimensional radiation attenuation prior tensor is mapped to the prior probability density function domain of the four-dimensional spatiotemporal interpolation grid;

[0044] The set of differential quotients of the dynamic photometric trajectory is converted into a likelihood observation vector; the joint posterior probability expectation deviation between the prior probability density function domain and the likelihood observation vector is solved to generate the spatiotemporal residual distribution matrix characterizing the degree of deviation between the two.

[0045] Further defining the step of mapping the three-dimensional radiation attenuation prior tensor to the prior probability density function domain of the four-dimensional spatiotemporal interpolation grid, the step further includes:

[0046] A virtual ray beam set is constructed with the camera's optical center as the origin; the in-path attenuation integral path of the virtual ray beam set as it penetrates the three-dimensional radiation attenuation prior tensor is calculated; based on the in-path attenuation integral path, expected values ​​and variance parameters are assigned to each time slice of the four-dimensional spatiotemporal interpolation grid to construct the prior probability density function domain.

[0047] Further defining step S4, the step of allocating fusion weight parameters and performing weighted fusion operation further includes:

[0048] Calculate the total number of valid samples contained in the set of differential quotients of the dynamic photometric trajectory; compare the total number of valid samples with the preset detection threshold.

[0049] When the total number of valid samples is greater than or equal to the preset detection threshold, a micro-dominant calculation strategy is executed: specifically, the state normalization parameters are extracted, mapped to a minimum weight value through a linear inverse proportional function, and assigned to the first sub-weight corresponding to the three-dimensional radiation attenuation prior tensor. At the same time, the difference between the value 1 and the first sub-weight is assigned to the second sub-weight corresponding to the set of differential quotients of the dynamic photometric trajectory. This ensures that when there are sufficient probes, the system is completely dominated by micro-observation data for judgment.

[0050] When the total number of valid samples is lower than the preset detection threshold, a downgrade calculation strategy is executed: the first sub-weight corresponding to the three-dimensional radiation attenuation prior tensor is monotonically increased, and the second sub-weight corresponding to the dynamic photometric trajectory differential quotient set is monotonically decreased, generating a positive detection indication signal based on the pure prior spectrum prediction mode.

[0051] Further defining the step of monotonically increasing the first sub-weight corresponding to the three-dimensional radiative attenuation prior tensor and simultaneously monotonically decreasing the second sub-weight corresponding to the set of differential quotients of dynamic photometric trajectories, the step further includes:

[0052] Obtain the total number of valid samples and the preset detection threshold.

[0053] Subtracting the preset detection threshold from the total number of valid samples yields a negative deviation value with a quantitative dimension.

[0054] The global residual norm of the spatiotemporal residual distribution matrix is ​​extracted, and a state normalization parameter is generated based on the global residual norm. It is hereby explicitly stated that the function of this parameter is to eliminate dimensional differences and to dynamically control the steepness of the weight degradation transition using the distribution state of the residuals.

[0055] The product of the negative deviation value and the state normalization parameter is calculated to obtain the dimensionless scaling deviation.

[0056] The scaling bias is input into a nonlinear relational mapping function based on data fitting (in this embodiment, a smooth variant of the Sigmoid function is used as a preset empirical algorithm model).

[0057] Further specifying the calculation logic for the nonlinear relational mapping function, the following steps are taken: using the natural constant as the base and the dimensionless scaling deviation as the exponent, perform an exponentiation operation to obtain an intermediate exponent result; add the value 1 to the intermediate exponent result to obtain the denominator; divide the value 1 by the denominator and extract the quotient as the first sub-weight; perform a subtraction operation to subtract the first sub-weight from the value 1 to obtain the second sub-weight.

[0058] Further defining step S5, the step of generating the reverse evolution indication signal further includes:

[0059] The spatiotemporal residual distribution matrix is ​​decomposed using a time-domain high-pass filter to isolate the high-frequency variation components with frequencies higher than the first cutoff frequency; the spatial clustering index of the high-frequency variation components is calculated.

[0060] When the spatial clustering index meets the preset convergence condition, all background noise nodes that do not meet the standard are blocked, and the average absolute value of the high-frequency variation component is obtained only for the high-frequency variation component that is marked as the core node.

[0061] A loss mapping transformation constant is introduced to eliminate the difference in physical dimensions, and the spatial clustering index, the absolute value of the average residual, and the loss mapping transformation constant are multiplied together to generate a dimensionless pure numerical result.

[0062] The dimensionless pure numerical result is established as the regularization penalty term, and the regularization penalty term is injected into the physical loss function of the graph neural network in an additive form. Then, the weight gradient of the hidden layer nodes of the network is calculated through the backpropagation algorithm, thereby forcing the graph neural network to update the local aerosol scattering coefficient, which is represented as the network node parameter, in the next calculation cycle.

[0063] Further specifying, the step of calculating the spatial clustering index of the high-frequency variation component further includes:

[0064] Determine the three-dimensional Euclidean distance matrix of the high-frequency variation component in the three-dimensional spatial coordinate system;

[0065] For each target node in the high-frequency mutation component, based on the three-dimensional Euclidean distance matrix, count the number of neighboring nodes whose distance is less than the preset search radius;

[0066] When the number of adjacent nodes is greater than a preset core density threshold, the target node is marked as a core node;

[0067] Summarize the total number of the marked core nodes;

[0068] By introducing a positive real number safety cushion parameter, the total number of nodes of the high-frequency variation component is added to the positive real number safety cushion parameter to obtain a safety denominator that is forced to be non-zero.

[0069] The ratio of the total number of core nodes to the safety denominator is calculated and used as the spatial clustering index.

[0070] The following is a detailed implementation description of the above content: An external structured data carrier, independent of the core computing logic, is established. This carrier is a configuration file in the device file system or a dataset table on a cloud server. Its technical function is to achieve hardware and software decoupling between the core algorithm logic and the application strategy. Operating parameters such as the astronomical ephemeris bias coefficient and the preset specular repulsion threshold are read from the external structured data carrier and extracted into the working memory of the current running process. The operating parameters are then assigned specific values. In this embodiment, the astronomical ephemeris bias coefficient is assigned a local actual logic compensation value based on the standard astronomical algorithm, which serves as the reference input constant for subsequent steps S1 to S3.

[0071] Astronomical ephemeris offset coefficient (denoted as ) This represents the logical compensation value of the solar altitude angle under latitude and longitude. Specifically, it reads and parses pre-approved geographical constants from a pre-built external structured data carrier.

[0072] In the description of the embodiments of the present invention, the "astronomical ephemeris bias coefficient" represents the logical compensation constant of the theoretical initial solar radiation energy extreme value under a specific three-dimensional spatial coordinate system, unaffected by meteorological disturbances in the lower troposphere; specifically, the logical derivation steps of the astronomical ephemeris bias coefficient are as follows:

[0073] Extract the absolute geographical latitude corresponding to the center point of the three-dimensional bounding box of the tunnel entrance; at the same time, extract the standard timestamp of the current calculation cycle, and based on the standard astronomical algorithm (the law of change of the obliquity of the ecliptic caused by the Earth's rotation and revolution), resolve the timestamp into the solar declination angle and local solar hour angle at the current moment.

[0074] Based on the spherical trigonometry theorem, the sine of the solar altitude angle at the current spatiotemporal node is calculated to quantify the absolute angle of incidence of sunlight relative to the ground plane. The specific calculation logic is as follows: the first intermediate component is obtained by multiplying the sine of the absolute geographical latitude and the sine of the solar declination angle; the second intermediate component is obtained by multiplying the cosine of the absolute geographical latitude, the cosine of the solar declination angle, and the cosine of the local solar hour angle; the first and second intermediate components are then added together, and the resulting sum is established as the sine of the solar altitude angle.

[0075] Since graph neural networks cannot directly process the solar constant with absolute physical dimensions, a network mapping reference constant is proactively introduced. The technical function of this reference constant is to proportionally map the absolute solar radiation constant at the top of the atmosphere into the normalized feature space of the graph neural network's hidden layers. The generation logic of the astronomical ephemeris bias coefficient is defined as follows: comparing the sine value of the solar altitude angle with zero, extracting the maximum value as the non-negative truncation value; multiplying the non-negative truncation value with the network mapping reference constant, the resulting product is established as the astronomical ephemeris bias coefficient. The core physical function of the above comparison operation logic for extracting the maximum value is: when the sun is below the horizon (nighttime condition, resulting in the calculated sine value of the solar altitude angle being less than zero), forcibly truncates the theoretical initial radiation energy extreme value to zero, to conform to the objective physical fact of passive optics.

[0076] The astronomical ephemeris bias coefficient, calculated through the above steps, plays a crucial technical role in the entire cross-scale assimilation calibration closed loop as an absolute upper limit for optical energy. This parameter is explicitly configured as the initial feature scalar of the "ephemeris radiation source node" in the macroscopic logical knowledge graph, providing physical boundary constraints for subsequent feature aggregation of network nodes simplified based on Maxwell's radiative transfer equations.

[0077] Preset specular repulsion threshold (denoted as) The threshold discreteness is used to determine the presence of invalid specular reflections on the target surface. The determination logic is as follows: A standard tunnel entrance is continuously acquired using an image acquisition interface with a resolution of at least 1080P and automatic gain control (AGC) disabled, under different sunlight angles. From the sample video streams, one thousand standard passenger vehicle sample sequences and one thousand special vehicle sample sequences with highly reflective metallic surfaces are extracted. The following coefficient of variation extraction operations are performed: For each desensitized dynamic target sample, all pixels within its bounding box are extracted, and the standard deviation and mean grayscale values ​​of the red, green, and blue color channels are calculated respectively. Since this scheme is applied to the dark area of ​​the tunnel entrance, low illumination may cause the mean grayscale value of the pixels to approach zero, thus causing the risk of division divergence. This embodiment introduces a preset dark field basis offset constant (set to a very small positive real number of 0.05). The technical function of this dark field basis offset constant is to ensure the non-zero lower limit of the denominator. The extracted mean grayscale values ​​of each channel are added to the dark field basis offset constant to obtain the corrected expected denominator of each channel. The standard deviation of the pixel grayscale value of each channel is divided by the corresponding corrected expected denominator to obtain the monochrome coefficient of variation of each channel. The arithmetic mean of the monochrome coefficients of variation of the three channels is performed to generate the absolutely convergent pixel coefficient of variation of the sample. A Gaussian distribution is fitted to the set of pixel variation coefficients of standard passenger vehicle samples; the 99th percentile value of this distribution curve is extracted; a preset safety margin constant of 0.05 is added to this percentile value to establish it as the preset specular rejection threshold. In this embodiment, the preset specular rejection threshold is assigned a value of 0.15. This calibration logic achieves a statistical balance between "preventing specular reflection from contaminating the photometric attenuation baseline" and "avoiding the mistaken elimination of normal bright light-colored vehicles," ensuring the physical consistency of the probe output in step S2.

[0078] Camera intrinsic transformation matrix (denoted as) The geometric distortion parameter, representing the mapping from two-dimensional pixel coordinates to three-dimensional world coordinates, is determined by extracting the factory calibration matrix of the monitoring equipment.

[0079] Preset detection threshold (denoted as) This represents the minimum number of samples required to maintain the confidence level of microscopic observations. Based on the law of large numbers in mathematical statistics, it is set as the minimum count value required to satisfy the 95% confidence interval. The logical derivation steps for the preset detection threshold in this embodiment are as follows: Obtain the variance of the optical probe extraction error of the target tunnel under historical benchmark operating conditions; simultaneously, based on the accuracy requirements for absolute brightness compensation, pre-calibrate the maximum tolerable relative error limit value. Based on the normal distribution interval estimation theory in mathematical statistics, extract the corresponding standard normal distribution critical value (a constant of 1.96 in this embodiment) for the set 95% confidence interval. Calculate the square of the standard normal distribution critical value and multiply it with the variance of the optical probe extraction error to obtain the intermediate numerator; calculate the square of the maximum relative error limit value as the denominator; divide the intermediate numerator by the denominator, round the calculated quotient up, and establish this rounded result as the preset detection threshold. In a preferred embodiment of the present invention, by substituting the calibration parameters of the tunnel entrance, the preset detection threshold is assigned to 500 probe individuals within 30 consecutive valid image frames, thereby constructing the basic environmental input constant benchmark of this embodiment.

[0080] First cutoff frequency (denoted as) The frequency used to characterize the time-domain boundary between slow environmental drift (regular gradual changes in sunlight and seasonal temperature variations) and sudden micrometeorological distortions (sudden dark fog and aerosol mutations caused by crosswinds) is determined by the following logical derivation steps: Extracting historical meteorological gradual change time-series data sequences from local meteorological monitoring stations within past continuous natural years, artificially marked as free of sudden severe weather. Applying a discrete-time to frequency-domain mathematical transformation algorithm (Fast Fourier Transform), converting the historical meteorological gradual change time-series data sequences into power spectral density sequences; performing discrete energy integration along the frequency axis starting from zero Hz on the power spectral density sequences until the accumulated integral value reaches 99% of the total energy of the sequence. The frequency coordinate point corresponding to the 99% total energy integration critical point is extracted and established as the first cutoff frequency. In a preferred embodiment of the invention, based on historical meteorological data statistics from valley terrain tunnels, this first cutoff frequency is explicitly calculated and assigned a value of 2.5 Hz. This serves as the absolute frequency division reference for a second-order infinite impulse response Butterworth time-domain high-pass filter, thereby isolating high-frequency variation components that characterize anomalous environmental changes.

[0081] The convergence rate of the prediction-observation residuals (denoted as ) The method represents the dynamic recovery capability after disturbance; it calculates the number of milliseconds required for the trace of the spatiotemporal residual distribution matrix to fall back to the steady-state baseline level after a sudden meteorological disturbance. Regarding the "prediction-observation residual convergence rate," the specific closed-loop calculation logic includes the following steps: within a historical, undisturbed operating cycle, extract the spatiotemporal residual distribution matrix for multiple consecutive calculation cycles, and calculate the trace of the matrix represented by the sum of the main diagonal elements of each matrix; perform an arithmetic mean operation on the values ​​of this series of traces to generate the expected value of the steady-state baseline; introduce a preset disturbance sensitivity factor, multiply the expected value of the steady-state baseline with this disturbance sensitivity factor, and generate the upper limit threshold for triggering dynamic disturbances. During real-time operation, the real-time trace value of the current spatiotemporal residual distribution matrix is ​​monitored synchronously; when the real-time trace value is greater than or equal to the upper limit threshold for triggering dynamic disturbances for the first time, a valid meteorological or photometric disturbance is determined, and the first system timestamp of the current moment is immediately extracted and locked. After locking the first system timestamp, the trace value of the spatiotemporal residual distribution matrix output in subsequent calculation cycles is continuously tracked. When the trace value shows a decaying trend and falls back to within the range of less than or equal to the steady-state baseline expectation plus the minimum tolerance constant for the first time, the disturbance is determined to have subsided, and the second system timestamp at the current moment is immediately extracted. The absolute value of the time difference between the second system timestamp and the first system timestamp is calculated to obtain the convergence duration in milliseconds. The reciprocal operation is performed, that is, the value is divided by the convergence duration, and the resulting quotient is established as the prediction-observation residual convergence rate. The disturbance sensitivity factor is limited to the real number range of (1.0, 5.0]. In this embodiment, the disturbance sensitivity factor is explicitly and preferably assigned a value of 2.0. The minimum tolerance constant is limited to the real number range of (0, 0.1). In this embodiment, the minimum tolerance constant is explicitly and preferably assigned a value of 0.05.

[0082] Cross-scale attention weight matrix (denoted as This represents the proportion of influence of macro-meteorology and micro-topography at different spatial nodes. The specific calculation and reasoning process is as follows:

[0083] Obtain a historical environmental dataset that meets the following key attributes: It includes a 10km-level meteorological cloud image tensor aligned to time resolution (as a global feature source) and a 100m-level 3D point cloud elevation model (as a local terrain source), and all input node data has undergone Z-score normalization to eliminate extreme value offsets. Perform cross-scale feature association calculations, specifically including the following steps: extract the corresponding dimension of the global meteorological feature map. The first feature vector, and the dimension corresponding to the local terrain structure map are... The second feature vector is introduced; a first feature transformation coefficient matrix and a second feature transformation coefficient matrix are introduced; these are used to map heterogeneous physical quantities to unified dimensionless values, so as to eliminate the difference in physical dimensions between meteorological attributes and spatial geometric attributes; matrix multiplication is performed on the first feature vector and the first feature transformation coefficient matrix to obtain a first target feature vector mapped to a unified dimensionless dimensional space; matrix multiplication is performed on the second feature vector and the second feature transformation coefficient matrix to obtain a second target feature vector mapped to the same dimensionless dimensional space; matrix dot product is performed on the first target feature vector and the second target feature vector to obtain an initial correlation score matrix; a normalization exponential function (Softmax function) is applied to each element in the initial correlation score matrix to compress the score to a preset normalization interval, thereby outputting the cross-scale attention weight matrix. The first feature transformation coefficient matrix is ​​defined as the mapping weight set from meteorological macroscopic attributes to the dimensionless latent space, and the data attribute is a two-dimensional real number matrix; in this embodiment, it is preferably obtained through pre-training. The first feature is a floating-point matrix; in practical applications, this parameter can take any tensor matching the feature dimension. The second feature transformation coefficient matrix is ​​defined as the mapping weight set from the terrain micro-geometric attributes to the dimensionless latent space, and the data attribute is a two-dimensional real matrix. In this embodiment, it is preferably obtained through pre-training. This is a floating-point matrix, and in practical applications, this parameter can take any tensor that matches the feature dimension.

[0084] Transient velocity scalar (denoted as) The velocity of the dynamic probe in three-dimensional space is determined by calculating the L2 distance of the motion optical flow vector between two consecutive frames. Specifically, the desensitized image sequence of the current frame and the adjacent previous frame is obtained. A dense optical flow estimation algorithm based on polynomial expansion is applied to extract the horizontal and vertical movement vectors of the central region of the same dynamic target bounding box, and the initial pixel displacement magnitude is obtained by performing a square root operation on them. The one-pixel diagonal length of the dynamic target bounding box in the previous frame and the current frame is extracted respectively. The pixel diagonal length of the previous frame is divided by the pixel diagonal length of the current frame to obtain a dimensionless relative depth compensation factor. The relative depth compensation factor uses the prior characteristic of constant vehicle physical size to characterize the depth scaling of perspective projection. The initial pixel displacement magnitude and the relative depth compensation factor are multiplied to perform spatial scale correction on two-dimensional perspective distortion, thereby outputting the transient velocity scalar with three-dimensional physical space consistency.

[0085] The attenuation integral path along the way (denoted as) This characterizes the cumulative loss of photons as they travel from the light source through regions of different aerosol densities to the viewpoint. It is determined by performing an accumulation operation on the optical thickness of each voxel node along the virtual ray beam direction. The logical derivation steps for the attenuation integral path along the path are as follows: obtain the retrieved camera intrinsic transformation matrix and extrinsic rotation and translation matrix; extract the absolute physical coordinates of the camera's three-dimensional optical center as the ray origin; for any reference target pixel in the current time slice of the four-dimensional spatiotemporal interpolation grid, calculate the world coordinate system direction vector corresponding to that pixel; combine the origin with this direction vector to establish a spatial virtual ray beam. Along the direction of the virtual ray beam, a fixed discrete step length parameter with absolute physical dimensions is set (0.5 meters in this embodiment). Starting from the origin, the sampling coordinate points are generated point by point outward according to the fixed discrete step length parameter, generating a series of sampling coordinate points evenly distributed in three-dimensional Euclidean space. For each sampling coordinate point, the nearest neighbor three-dimensional voxel node is searched in the three-dimensional radiation attenuation prior tensor generated in step S1, and the prior attenuation feature scalar stored in the voxel node is extracted. Here, the prior attenuation feature scalar is explicitly defined as the tensor element value of the corresponding spatial position derived by the graph neural network in step S1 based on the aggregation of local aerosol scattering coefficient and local absorption attenuation coefficient. The prior attenuation feature scalar corresponding to any of the above-obtained sampling coordinate points is extracted; the prior attenuation feature scalar is multiplied by the fixed discrete step length parameter to calculate and generate the discrete optical thickness within the single step length region; an algebraic accumulation operation is performed on all discrete optical thicknesses generated along the path of the virtual ray beam. The final sum output by the above traversal and accumulation operations is physically equivalent to the total cumulative optical loss of the photon propagating along the ray. The corresponding spatial mapping index chain generated by this calculation process is then established as the attenuation integral path along the way.

[0086] Core density threshold (denoted as) The lower limit number of nodes is used to determine whether local high-frequency variations constitute a clustered meteorological anomaly. The determination logic is as follows: it is set based on the product of the total grid volume within the preset search radius and the empirical fog density. Specifically, the preset search radius parameter is obtained. In this embodiment, for the physical motion limit of the matching vehicle at the 2.5 Hz cutoff frequency, the preset search radius parameter is assigned to the actual physical distance (preferably 30 to 50 meters) corresponding to the step size of three consecutive nodes in the three-dimensional tetrahedral interpolation grid. With the target node as the center, based on the three-dimensional physical side length of a single voxel, the total number of discrete grid nodes completely contained in the three-dimensional Euclidean space enveloped by the search radius parameter is counted, and it is forcibly mapped to the dimensionless local search space volume. The spatiotemporal residual distribution matrix is ​​obtained under ideal clear sky and interference-free historical operating conditions with a continuous operating time of not less than 72 hours. Under these ideal operating conditions, the discrete time average of the high-frequency variation node density caused by the equipment's background thermal noise per unit volume within the 72-hour time window is calculated and established as the mathematical expectation value of the density. The calculated local search space volume is multiplied by the expected value of the background thermal noise density to obtain the number of baseline nodes. A significance confidence multiplier of 3.0 is introduced, and the number of baseline nodes is multiplied by the significance confidence multiplier. The result is then rounded down to generate the final core density threshold. This ensures that only when the spatial aggregation of high-frequency variations reaches a statistically significant level, representing a real local exhaust gas agglomeration or aerosol mutation, will a node be identified as a core node, thus eliminating the possibility of erroneous reverse evolution triggered by random white noise in the environment.

[0087] Further explanation of steps S1 to S5 above:

[0088] For step S1, in this embodiment, the "tunnel entrance transition space" is a three-dimensional topological region with physical constraints and algorithmic boundaries. Specifically, the physical boundary of this space is uniquely defined by the following three dimensions of features:

[0089] The first dimension is the optical and structural physical boundary: this space encompasses a three-dimensional physical region extending from the end of the tunnel's external approach section to the tunnel's internal entrance section and brightness transition section. Essentially, it is a geometric set representing the nonlinear, dramatic change in natural ambient illuminance to artificial lighting illuminance—the actual physical domain where the "black hole effect" occurs.

[0090] The second dimension is the hardware-sensing field of view boundary: the boundary of this space on the two-dimensional projection plane, constrained by the effective monitoring field of view projection cone jointly determined by the "preset camera intrinsic parameter transformation matrix" and extrinsic parameter matrix in step S3. Any physical region exceeding the penetration range of this set of virtual ray beams is not considered a processing object of this method.

[0091] The third dimension is the algorithm grid topology boundary: at the algorithmic logic level, this space is discretized and equivalent to the "four-dimensional spatiotemporal interpolation grid" constructed in step S3. The limit effective distance of its depth coordinate axis depends on the maximum spatial displacement vector traversed by the dynamic target in step S2 when satisfying the "preset effective tracking frame length".

[0092] In this embodiment, the "tunnel entrance transition space" is a deterministic computational domain formed by the intersection of macroscopic geographical topography, mesoscopic camera grating projection field of view, and microscopic aerodynamic attenuation law. The absolute physical coordinate range of this space is determined based on the on-site camera calibration parameters and standard tunnel lighting design drawings.

[0093] The specific construction and instantiation process for the step of "acquiring macroscopic logical knowledge graph data containing target geographic coordinates, ephemeris status, and real-time meteorological cloud cover status" further includes the following operations:

[0094] The absolute geographic coordinates (longitude, latitude, and elevation) of the center point of the pre-defined 3D spatial bounding box at the tunnel entrance are extracted. Based on these absolute geographic coordinates, real-time meteorological cloud cover tensors (including space meteorological elements such as cloud thickness and aerosol optical thickness) and local micro-topographic point cloud data at the corresponding timestamps are synchronously acquired through a standardized application programming interface (API). These heterogeneous data are abstracted into two types of core entity nodes in the atlas data structure. One is the "meteorological environment node," which maps the discrete grid data in the real-time meteorological cloud cover tensor to the hidden layer attribute vector of this type of node; the other is the "topographic boundary node," which maps the local micro-topographic point cloud data to attribute vectors representing the geometric shape of spatial obstructions. The astronomical ephemeris bias coefficients retrieved from the working memory are explicitly configured as "ephemeris radiation source nodes" with globally broadcast attributes, whose characteristic scalars strictly characterize the theoretical limit of incident solar radiation intensity at the current time and latitude / longitude.

[0095] Relationship edges are constructed using optical physical spatial proximity, specifically based on the photon rectilinear propagation laws and spatial adjacency matrix of the real three-dimensional world. Directed connections are established between "ephemeral radiation source nodes" and "meteorological environment nodes," as well as between "meteorological environment nodes" and "terrain boundary nodes." Each directed connection naturally represents the energy transfer medium channel through which photons undergo absorption or scattering physical interactions along that spatial path. The dataset containing ephemeral radiation source nodes, meteorological environment nodes, terrain boundary nodes, and physically directed connections is collectively instantiated into macroscopic logical knowledge graph data.

[0096] Furthermore, this process retrieves astronomical ephemeris bias coefficients from memory. It extracts macroscopic logical knowledge graph data containing the target geographic coordinates ephemeris status and real-time meteorological cloud cover status. A graph embedding algorithm is used to perform spatial rasterization mapping on entity nodes and relational edges in the macroscopic logical knowledge graph data to generate a global meteorological feature map and a local terrain structure map with a uniform grid resolution. Graph neural network topology inference is executed to address the scale discontinuity between 100-kilometer-level meteorological data and 100-meter-level terrain data. Specifically, a graph neural network topology inference process containing at least three layers of graph convolutional network (GCN) layers is executed, with the GCN layers using the ReLU activation function. Based on Maxwell's radiative transfer equation, the simplified boundary condition constraints are injected as a regularization penalty term into the loss function of the graph neural network for iterative training to calculate the three-dimensional radiative attenuation prior tensor with spatial depth attributes. Specifically, the global meteorological feature map and local terrain structure map are extracted from the macroscopic logical knowledge graph data; a cross-scale attention weight matrix between them is calculated. Using this cross-scale attention weight matrix, a weighted aggregation operation is performed on the initial three-dimensional radiation field. Then, the physical constraint logic based on the simplification of Maxwell's radiative transfer equation is executed: specifically for any hidden layer node in the graph neural network, the received incident radiation feature scalar and the output outgoing radiation feature scalar are extracted.

[0097] Regarding the "local aerosol scattering coefficient," to avoid the randomness and blindness of network weights in the cold start state, the initial value generation logic follows these steps before the initial computation cycle of the graph neural network topology inference is executed for the first time:

[0098] From the acquired real-time meteorological cloud cover data, the initial visibility scalar and particulate matter concentration scalar corresponding to the current target spatial grid node are extracted. A pre-calibrated atmospheric extinction empirical transformation constant (a constant parameter derived based on the classical Koschmieder law) is introduced; the baseline extinction coefficient is obtained by dividing 3.912 times the natural logarithm constant by the initial visibility scalar; the baseline extinction coefficient is then linearly scaled and weighted using the particulate matter concentration scalar to calculate the initial physical scattering baseline value that purely reflects the physical atmospheric properties. The calculated initial physical scattering baseline value is then explicitly assigned to the "learnable spatial weight parameter" corresponding to a specific hidden layer in the graph neural network. The local aerosol scattering coefficient and local absorption attenuation coefficient corresponding to the spatial location of the node are extracted. A physical mapping transformation constant with a value range limited to (0,1) is actively introduced. In this embodiment, the preferred value of the physical mapping transformation constant is 0.65. The role of this physical mapping transformation constant is to eliminate the dimensional difference between the abstract features of the neural network and the absolute optical energy. The incident radiation feature scalar is multiplied by the physical mapping transformation constant to obtain a standardized incident energy value. This standardized incident energy value is then multiplied by the local aerosol scattering coefficient and the local absorption attenuation coefficient respectively, and the sum is obtained to obtain the theoretical total loss value. The absolute difference between the standardized incident energy value minus the theoretical total loss value and the emitted radiation feature scalar is calculated. This absolute difference of all nodes is accumulated as a physical regularization penalty term and injected into the network, thereby outputting a three-dimensional radiation attenuation prior tensor that follows physical conservation and has spatial depth attributes. This mechanism ensures that the microscopic truncation effect of complex terrain on large-scale meteorological illumination is accurately quantified.

[0099] For step S2, after obtaining the original image sequence in this embodiment, an irreversible edge anonymization and desensitization operation is performed by applying a kernel size not less than [amount missing]. And the standard deviation parameter A Gaussian blur mask of at least 3.0 is used to perform irreversible pixel obfuscation on the target area to remove license plate and face identifiers. Dynamic target bounding box sequences are extracted from the desensitized images. Specifically, dynamic target bounding box sequences satisfying a preset effective tracking frame length are extracted from consecutive desensitized temporal images. The preset effective tracking frame length represents the minimum number of effective frames for continuous and stable tracking of the dynamic target in the field of view (preferably at least 15 frames in this embodiment).

[0100] Calculate the monochrome variation coefficients across multiple color channels in the first color space, and perform an arithmetic mean operation on these monochrome variation coefficients to generate the pixel variation coefficients of the dynamic target. Compare these pixel variation coefficients with a preset specular repulsion threshold. For valid probes with a value below the preset specular repulsion threshold, calculate their apparent grayscale temporal difference.

[0101] If a vehicle passes through a tunnel entrance at 10 km / h or 100 km / h, the grayscale differences within the same time window lack physical comparability. Therefore, dimensional normalization is performed, specifically extracting the motion optical flow vectors of consecutive target frames, and calculating the transient velocity scalar. A dynamic normalization parameter is introduced, and the transient velocity scalar is multiplied by the dynamic normalization parameter to generate a dimensionless velocity penalty factor. The ratio of the apparent gray-scale temporal difference to the velocity penalty factor is calculated, and a mapping function characterized by the Sigmoid nonlinear activation function is introduced to map the ratio to a dimensionless preset normalization interval of [0,1], thereby generating a set of dynamic photometric trajectory differential quotients characterizing the optical properties of the environment.

[0102] Regarding the process of "calculating the ratio of apparent gray-scale temporal difference to transient velocity scalar," this embodiment specifically implements it through the following sub-steps: Obtaining the extracted apparent gray-scale temporal difference and transient velocity scalar; since the apparent gray-scale temporal difference and transient velocity scalar belong to completely different physical properties, a dynamic normalization parameter is introduced; this makes the photometric temporal difference at different motion velocities comparable; multiplying the transient velocity scalar with the dynamic normalization parameter to obtain a dimensionless velocity penalty factor; dividing the apparent gray-scale temporal difference by the velocity penalty factor to obtain the ratio. In this embodiment, the dynamic normalization parameter is defined as a scaling factor to eliminate the heterogeneous dimensions of velocity and gray-scale, and is preferably... In practical applications, this parameter takes any real number between (0, 0.1).

[0103] In another embodiment, for obtaining the apparent grayscale temporal difference, the average grayscale value of the dynamic target bounding box in the current frame and the average grayscale value in the previous frame are obtained; a subtraction operation is performed on the two and the absolute value is taken to obtain the apparent grayscale temporal difference. A preset reference speed value is obtained; this reference speed value is a kinematic scalar characterizing the expected kinematic speed of a vehicle passing through the tunnel entrance under ideal free-flow conditions. Its determination logic is as follows: the legally approved maximum design speed of the target tunnel entrance section under traffic management regulations is extracted, and based on the camera's calibration frame rate and spatial physical resolution, this legally approved maximum design speed is mathematically converted and uniformly mapped to a pixel-level reference scalar with the same dimension as the transient speed scalar.

[0104] The ratio of the acquired transient velocity scalar to the reference velocity value is calculated to obtain the relative velocity ratio. A mapping weighting factor with a value range of (0,1) is introduced, which maps heterogeneous physical quantities to uniform dimensionless values ​​to achieve parameter normalization. The relative velocity ratio is multiplied by the mapping weighting factor to generate a velocity attenuation penalty factor, thereby obtaining the apparent gray-scale temporal difference of the same dynamic target on continuous time slices.

[0105] Introducing a minimal positive real constant The speed attenuation penalty factor is added to this minimal positive real constant to obtain a corrected divisor term. The apparent gray-scale temporal difference is then divided by the corrected divisor term to eliminate the interference of target motion speed differences on calibration, thereby outputting a set of dynamic photometric trajectory differential quotient elements with equivalent comparability. This ensures that the photometric attenuation data output by vehicles acting as "optical probes," whether crawling at low speeds or entering tunnels at high speeds, has absolute equivalence and comparability in the spatiotemporal topological space. This step fundamentally transforms vehicles traveling on the road into "dynamic optical probes," eliminating the need for external standard light sources and professional luminance meters. An orthogonal verification mechanism for pixel variation coefficients is introduced to mathematically eliminate baseline drift contamination caused by complex reflective materials. This makes the anti-interference capability of microscopic photometric sampling stronger the traffic flow.

[0106] For step S3, based on the retrieved camera intrinsic transformation matrix and extrinsic rotation and translation matrix, and combined with a preset time sliding window, a four-dimensional spatiotemporal interpolation grid covering the real physical coordinates and time dimension of the three-dimensional world is constructed. A set of virtual ray beams with the camera's optical center as the origin is constructed to solve the projection distortion of the three-dimensional spatial field onto the 2D camera's viewpoint. The in-path attenuation integral path of this ray beam as it penetrates the three-dimensional radiation attenuation prior tensor is calculated. Based on the accumulated optical loss along this path, the expected value and variance are assigned to the time slices of the four-dimensional grid, and the prior probability density function domain is constructed. In a preferred embodiment, the process of "assigning expected value and variance parameters to each time slice of the four-dimensional spatiotemporal interpolation grid based on the in-path attenuation integral path" is specifically implemented through the following sub-steps:

[0107] The cumulative optical loss along the attenuation integral path is obtained and directly established as the expected value of a normal distribution for each time slice. The weather forecast confidence interval attached to the real-time meteorological cloud cover data obtained in the previous step is extracted to determine the uncertainty measure of the probability distribution. An uncertainty weighting factor is introduced; this factor maps the dimensionless confidence index of the weather forecast to an absolute uncertainty value with optical attenuation dimensions, achieving a smooth transition between heterogeneous parameter domains. The median of the weather forecast confidence interval is subtracted from the numerical value 1 to obtain the weather fluctuation error rate. The weather fluctuation error rate and the uncertainty weighting factor are multiplied to calculate the discrete standard deviation parameter. The discrete standard deviation parameter is squared to generate the variance parameter. Based on the established expected value and variance parameter, a prior probability density function domain is generated in a four-dimensional spatiotemporal interpolation grid.

[0108] The median of the weather forecast confidence interval is defined as a characterization of the accuracy of external macro-meteorological data sources and is a dimensionless percentage constant; in this embodiment, it is preferably 0.85, and in practical applications, this parameter can be any real number between [0.5, 0.99]. The uncertainty weighting factor is defined as the mapping multiplier from the meteorological error rate to the optical variance; in this embodiment, it is preferably 10.0, and in practical applications, this parameter can be any real number between [1.0, 50.0].

[0109] The dimensionless set of differential quotients is converted into likelihood observation vectors. Specifically, based on Bayes' theorem, multiplication and normalization constant division operations are performed. By calculating the Koolbek-Leibler divergence (KL divergence), the joint posterior probability expectation deviation between the two is solved, and the spatiotemporal residual distribution matrix of the mapping error is output. The dynamic performance of this process is characterized in real time by the convergence rate of the prediction-observation residuals. The steps in this embodiment to convert the set of dynamic photometric trajectory differential quotients into likelihood observation vectors specifically include the following spatial registration mapping operations:

[0110] The occurrence timestamps and three-dimensional spatial physical coordinates of each discrete element in the set of differential quotients of dynamic photometric trajectories are extracted. For each discrete element, at least three adjacent grid reference nodes with the closest Euclidean distance are searched in a pre-constructed four-dimensional spatiotemporal interpolation grid. A dimensionless spatial distance attenuation factor is introduced, which is calculated based on the absolute physical distance between the discrete probe and the grid node. The differential quotient values ​​of the discrete elements are multiplied by the spatial distance attenuation factors corresponding to each adjacent grid reference node to obtain multiple projection components. All projection components received by each reference node in the four-dimensional spatiotemporal interpolation grid are accumulated and local mean smoothed to transform the disordered discrete trajectory lattice into a likelihood observation vector with the same dimensional structure as the prior tensor.

[0111] For step S4, obtain the set of differential quotients of the input dynamic photometric trajectory; for all discrete elements in this set, perform the operation based on the Laida criterion (3). The statistical outlier removal operation of the criteria is used to physically isolate singular differential quotient elements caused by occasional blindness due to strong light or abnormal shaking of the probe vehicle; the number of element individuals that are actually left after the removal operation is counted and established as the total number of effective samples contained in the dynamic photometric trajectory differential quotient set.

[0112] The result is compared with a preset detection threshold. If a downgrade condition is triggered, the difference between the two is calculated. This difference is then input into a smooth inverse proportional activation function (a variant of the sigmoid function). The output of this activation function is limited to the interval (0.5, 1.0). This output value is extracted as the first sub-weight representing the prior weight; a subtraction operation is performed, subtracting the first sub-weight from the value 1 to obtain the second sub-weight representing the observation weight. Based on these two normalized weights, a weighted addition operation is performed to generate a positive detection indication signal.

[0113] This embodiment, based on fusion weight parameters, performs element-wise tensor weighted addition operations on the three-dimensional radiation attenuation prior tensor (as the first matrix) and the set of dynamic photometric trajectory differential quotients transformed into isomorphic likelihood observation vectors (as the second matrix). This generates a three-dimensional brightness attenuation compensation control matrix with precise spatial coordinate mapping properties, and encapsulates it into a positive detection indication signal characterizing the global adaptive calibration absolute brightness attenuation distribution. This signal is then sent to the tunnel lighting programmable logic controller via the industrial control bus for segmented closed-loop dimming control.

[0114] Specifically, the total number of valid samples and the preset detection threshold value are obtained; the total number of valid samples is subtracted from the preset detection threshold value to obtain the negative deviation value (this value is negative because the degradation condition requires the total number of samples to be lower than the threshold value).

[0115] Perform inverse proportional nonlinear activation operation: For the step of extracting state normalization parameters, obtain the values ​​of all elements in the spatiotemporal residual distribution matrix; perform Frobenius norm calculation logic on the spatiotemporal residual distribution matrix, specifically calculating the sum of squares of all element values, and performing a square root operation on the sum of squares to extract the global residual norm of a single scalar; introduce a preset degradation sensitivity adjustment constant; multiply the global residual norm and the degradation sensitivity adjustment constant, and the resulting product is established as the state normalization parameter. Multiply the negative bias value with the calculated state normalization parameter to obtain the scaling bias;

[0116] Using the natural constant as the base and the scaling deviation as the exponent, a power operation is performed to obtain the exponent result. The value 1 is added to the exponent result to obtain the denominator. The value 1 is then divided by the denominator. The quotient of this division operation is extracted. Since the scaling deviation is negative, the quotient mathematically converges to the range of 0.5 to 1.0. This quotient is directly established as the first sub-weight. This ensures smooth degradation through the intrinsic properties of the mathematical function, eliminating nonlinear oscillations or even equipment damage to the actuator (tunnel lighting controller) caused by step signals.

[0117] The degradation sensitivity adjustment constant is defined as a tolerance adjustment factor for the distribution of meteorological forecast errors. This constant has a joint physical dimension that is the reciprocal of the product of the negative deviation value and the global residual norm. Its core technical function is to cancel out the combined dimension of the sample size and optical deviation in the multiplication operation, ensuring that the output scaling deviation is a purely dimensionless scalar. In this embodiment, its value is preferably 0.1. In practical applications, this parameter can take any real number between (0, 1.0).

[0118] This step focuses on the extremely high robustness of the forward detection loop. By explicitly performing Frobenius norm calculation and introducing a degradation sensitivity adjustment constant, feature loss during matrix-to-scalar dimensionality reduction is avoided, enabling the precise quantization of multi-dimensional spatiotemporal residuals into a single control variable. When the observation end experiences physical failure due to force majeure (obscuration of the horizon or loss of the observation object), the entire computational chain smoothly transitions to a "pure prior map prediction mode" due to the existence of the aforementioned nonlinear fusion weight parameter. This improves the stability of strategy triggering.

[0119] For step S5, a second-order infinite impulse response (IIR) Butterworth time-domain high-pass filter is applied to decompose the spatiotemporal residual distribution matrix. Utilizing its characteristic of having the flattest amplitude-frequency response within its passband to avoid introducing phase distortion, specifically, along the independent time axis dimension of the spatiotemporal residual distribution matrix (four-dimensional tensor), a one-dimensional signal array sliding window extraction and frequency decomposition are performed on each fixed three-dimensional physical node to isolate high-frequency variant components with frequencies higher than the first cutoff frequency. The three-dimensional Euclidean distance matrix of this component in the three-dimensional spatial coordinate system is calculated to eliminate interference from random environmental white noise. The number of core nodes with node density exceeding the core density threshold within the preset search radius is counted. A division operation is performed to prevent the hardware divergence and collapse of the division arithmetic logic unit (ALU) caused by the total number of high-frequency variant components approaching zero due to extreme global environmental stability. A positive real number safety cushion parameter (assigned to a value in this embodiment) is introduced. ); perform an addition operation on the total number of nodes of the high-frequency variation component and the positive real number safety cushion parameter to obtain a forced non-zero safety denominator; calculate the ratio of the total number of core nodes to the safety denominator as the spatial clustering index characterizing the absolute convergence state.

[0120] For determining whether the spatial clustering index meets the preset convergence condition: read the preset critical clustering threshold from the external structured data carrier; determine whether the current spatial clustering index is greater than or equal to the critical clustering threshold.

[0121] If the judgment result is negative, the error distribution is determined to be pure random white noise, the injection penalty term is abandoned, and the next regular calculation cycle begins.

[0122] If the judgment result is yes, then the current error feature is determined to meet the preset convergence condition (exhibiting non-random physical clustering features, including sudden dark fog), and the adaptive feedback adjustment mechanism is immediately started to obtain the average absolute value of the high-frequency variation component, and simultaneously obtain the convergence rate of the prediction-observation residual generated by the previous monitoring.

[0123] A pre-calibrated loss mapping transformation constant is introduced, with its value limited to the interval (0,1] (preferably 0.01 in this embodiment). The technical function of this loss mapping transformation constant is defined as eliminating the difference in physical dimensions between spatial clustering density, optical residuals, and neural network optimization gradients. Then, dynamic feedback fusion calculation logic is executed: the spatial clustering index, the absolute value of the average residual, and the loss mapping transformation constant are multiplied together to generate an initial penalty benchmark value; the prediction-observation residual convergence rate is introduced as the feedback denominator. To avoid the risk of division by zero divergence during the calculation process, a positive real number safety cushion parameter is introduced, and the prediction-observation residual convergence rate is multiplied by this positive real number. The safety cushion parameters are added together to obtain the corrected convergence rate. The initial penalty benchmark value is divided by the corrected convergence rate to generate a dimensionless pure numerical result that is dynamically amplified, and this pure numerical result is established as the final regularization penalty term. The regularization penalty term is injected into the physical loss function of the graph neural network in additive form, and then the weight gradient of the hidden layer nodes of the network is calculated through the error backpropagation algorithm. The lower the convergence rate of the perturbed residual, the higher the dynamic amplification factor of the calculated output regularization penalty term, thereby forcing the graph neural network to update the local aerosol scattering coefficient with a larger step size, achieving a highly sensitive self-healing loop closure across evolutionary cycles.

[0124] Before the next computation cycle begins, the regularization penalty term is injected additively into the physical loss function of the graph neural network, and the specific topology of the graph neural network is forced to be declared as a Physics-Informed-Graph-Neural-Network. This explicitly binds the learnable spatial weight parameters representing the attenuation ratio of information transmission between nodes in a specific hidden layer to the corresponding physical-level local aerosol scattering coefficients. Then, after calculating the adjustment gradient for updating these learnable spatial weight parameters using the gradient descent optimizer, a network parameter clamping operation based on the physical conservation boundary is enforced, specifically including:

[0125] After calculating the unconstrained candidate update weights, a physical boundary activation function (ReLU or a limiting function with a lower bound) with non-negative interval mapping characteristics is connected in series at the network output of the hidden layer. The unconstrained candidate update weights are extracted and input into the physical boundary activation function. If the calculated update gradient is too large due to a severe regularization penalty, causing the candidate weight to exhibit a negative value that violates objective physical laws, the physical boundary activation function forcibly truncates it to a minimal non-negative positive real number (10 to the power of -4 in this embodiment). The output value after non-negative activation constraint is equivalently established as the local aerosol scattering coefficient in the next calculation cycle, thus outputting a reverse evolution indicator signal. This mechanism ensures that regardless of the dynamically amplified penalty gradient caused by the low residual convergence rate during the reverse evolution calculation cycle, the final output local aerosol scattering coefficient is always absolutely constrained within the true non-negative natural meteorological physical extreme boundary.

[0126] The critical clustering threshold is defined as the statistical confidence threshold for distinguishing random environmental thermal noise from real physical meteorological clustering anomalies. In this embodiment, it is preferably 0.15, but in practical applications, this parameter can be any real number between [0.10, 0.30]. The fuzzy convergence condition is visualized as a comparison logic with the critical clustering threshold; the above is further elaborated as follows:

[0127] Using a nonlinear relational mapping function based on data fitting, the system outputs a first sub-weight assigned to the three-dimensional radiative attenuation prior tensor and a second sub-weight assigned to the set of differential quotients of dynamic photometric trajectories. Following the normalization principle, a forced subtraction operation is performed, subtracting the first sub-weight from the value 1 to obtain the second sub-weight, ensuring that the sum of the two is always equal to the value 1.

[0128] When the first sub-weight approaches its theoretical minimum (and the second sub-weight approaches 1), the representation is in an "absolutely microscopic dominant state." In this physical state, the tunnel entrance area has traffic flow and is not affected by severe weather-induced blindness; the number of effective probes within the set of dynamic photometric trajectory differential quotients fully meets the requirements of the law of large numbers in statistics. Adaptive calibration is performed entirely based on high-frequency, high-resolution real-time dynamic photometric trajectories to output the highest level of spatial brightness compensation control accuracy. When the first sub-weight approaches its theoretical maximum (converging to 1 in this embodiment), the representation encounters an "observation source cliff-like failure state." This state corresponds to fog, heavy rain obscuring the view, or zero traffic flow at night, where the number of microscopic photometric probes is below the preset detection threshold required to maintain the confidence level of microscopic observations. Under this trend, a smooth and forced switch to a "pure prior map prediction mode" is implemented, relying entirely on a three-dimensional radiation attenuation prior tensor containing historical macroscopic meteorological and topographic features to generate a positive detection indication signal. This trend's convergent design avoids the phenomenon of command divergence or system crashes in industrial control systems caused by the loss of input information entropy, and establishes the absolute safety of the physical boundaries of the lighting control system.

[0129] Furthermore, regarding the correlation analysis between the total effective sample size and the first sub-weight: within the interval triggering the degradation calculation strategy, the total effective sample size and the first sub-weight exhibit a non-linear negative correlation. When the number of available sensors or probes decreases linearly, the statistical error caused by insufficient sampling is non-linearly amplified. A negative bias value is generated by calculating the difference between the total effective sample size and the preset detection threshold, and the scaling bias containing this negative bias value is used as the exponent of the natural constant for exponential calculation. This negative mapping structure precisely matches the objective law of uncertainty growth when information sources are depleted. This negative correlation design ensures that when the total effective sample size just exceeds the preset detection threshold, the weight transition is in a smooth evolution period, physically avoiding control command jumps caused by occasional traffic flow fluctuations; while when effective probes are extremely scarce, the weight rapidly tilts towards the prior end, supporting the beneficial effect of "providing an absolutely smooth safety control degradation curve under extreme failure conditions".

[0130] Correlation analysis of global residual norm and weight degradation steepness: The global residual norm extracted from the spatiotemporal residual distribution matrix is ​​positively correlated with the gradient steepness of the nonlinear relational mapping function based on data fitting. The global residual norm quantitatively characterizes the overall physical deviation between macroscopic prior predictions and microscopic probe observations. The larger the norm value, the more severe the deviation of the current micro-meteorological environment (local aerosol distortion caused by sudden crosswinds) from the macro-weather forecast. Under this objective premise, if the number of microscopic probes falls below a threshold, it reverts to the prior state at a faster convergence rate, physically severing the severely contaminated erroneous visual observation data link. Introducing the global residual norm into the calculation process of the state normalization parameter enables the input variables of the nonlinear activation function to have environmental adaptive adjustment capabilities. This operational logic completely abandons the rigid fixed degradation curve, realizing a dynamic optimization mechanism of "faster degradation for larger residuals and slower degradation for smaller residuals," maximizing control robustness under complex weather and traffic flow disturbances.

[0131] Furthermore, this embodiment is configured in the following digital twin verification monitoring scenario: a tunnel entrance environment located on a mountainous highway. This physical scenario is subject to dual independent disturbances from high-frequency micro-meteorological changes (random generation and dissipation of valley advection fog) and irregular traffic flow attenuation. The specific parameter settings of the mathematical model are as follows:

[0132] The preset detection threshold is set at 500 discrete individuals; the degradation sensitivity adjustment constant is fixed at 0.1; the three-dimensional radiation attenuation prior tensor output by macroscopic logic knowledge graph reasoning maintains a stable output within the current verification cycle.

[0133] This table aims to verify the technical advantages of the dynamic weight allocation mechanism of this embodiment compared to existing technologies when encountering varying degrees of failure of microscopic information sources. Existing tunnel brightness compensation technologies generally employ rigid Boolean logic judgments, that is, when the number of effective probes falls below a preset limit, the observation source is immediately cut off, causing a precipitous jump in control weights. This table introduces a "traditional Boolean hard switching scheme" under the same operating conditions as a comparative reference, recording the conversion process of the input features and internal scaling deviation of this invention into the normalized first sub-weight, quantitatively revealing the beneficial effect of this invention in absorbing the impact of sudden environmental changes and achieving smooth control degradation by relying on the intrinsic properties of mathematical functions.

[0134] Table 1: Calculation Examples of System Response Comparison under Different Environmental Disturbances

[0135] Scene state test points Total number of valid samples Global residual norm Scaling deviation First weight of traditional Boolean scheme First sub-weight of the present invention Step impact reduction Regular traffic flow and stable weather 800 1.5 Downgrade not triggered 0.05 0.05 0 Critical traffic flow and weak disturbance 490 1.2 -1.2 1 0.768 0.232 Traffic flow decreased and moderate disturbance 470 3 -9 1 0.999 0.001 Severe congestion and severe distortion 450 6.5 -32.5 1 1 0 Extremely low flow and severe distortion 400 6.5 -65 1 1 0 Pure nighttime empty period with zero traffic 50 2 -90 1 1 0

[0136] The derivation logic for setting the step impact reduction magnitude is as follows: obtain the first weight of the traditional Boolean scheme at the current test point; obtain the first sub-weight of the present invention at the current test point; subtract the first sub-weight of the present invention from the first weight of the traditional Boolean scheme, and extract the absolute difference as the step impact reduction magnitude.

[0137] The step impact reduction magnitude quantitatively characterizes the amount of transient control command deviation absorbed and offset by the nonlinear activation logic of this invention compared to the traditional threshold cutoff strategy within the critical range triggering degradation. A larger value indicates a more significant hardware protection effect in avoiding mechanical oscillations in industrial actuators.

[0138] Based on the data set shown in Table 1, a representative critical failure scenario of "critical traffic flow with weak disturbance" was extracted for in-depth technical analysis. In this physical scenario, the number of effective samples (490) falling below the preset detection threshold (500) is relatively small. According to the traditional Boolean scheme judgment logic of existing technology, since the total number of effective samples triggers the threshold, the highest level of degradation response is forcibly initiated, and its first weight self-stabilizing state of 0.050 experiences a precipitous increase, instantly jumping to the extreme value of 1.000. This absolute step deviation of 0.950 will be directly converted into a destructive step voltage or pulse control signal in actual industrial control buses, directly inducing nonlinear oscillations in the tunnel lighting programmable logic controller.

[0139] Conversely, this invention performs multi-source parameter smoothing fusion control. The total effective sample size in the scenario is extracted and subtracted from a preset detection threshold, yielding a negative bias value of -10 with a quantitative dimension. Simultaneously, the global residual norm (1.20) extracted from the four-dimensional spatiotemporal interpolation grid mapping error is multiplied by a degradation sensitivity adjustment constant with a fixed value of 0.1, generating a state normalization parameter of 0.12. The negative bias value is then multiplied by the state normalization parameter to calculate the scaling bias that completely eliminates the composite physical dimensions, which has a value of -1.20.

[0140] Entering the nonlinear mapping stage, the natural constant is used as the base, and the calculated -1.20 is used as the exponent to perform a power operation, generating an intermediate exponent result. The value 1 is added to this intermediate exponent result to obtain the denominator, and finally, the value 1 is divided by this denominator. Through this smooth mapping based on the intrinsic properties of the function, the output first sub-weight of this invention naturally converges to 0.768, without exhibiting extreme binary jumps.

[0141] The absolute jump variable of the traditional scheme when the threshold is broken is extracted (the difference between 1.000 and 0.050 equals 0.950); then the absolute jump variable of the present invention under the same conditions is extracted (the difference between 0.768 and 0.050 equals 0.718). The difference between the absolute jump variable of the traditional scheme and the absolute jump variable of the present invention (0.232) is obtained, and this difference is divided by the absolute jump variable of the traditional scheme. Under critical micro-meteorological disturbance conditions, the present invention reduces the transient step impact amplitude of the control command by 24.4%.

[0142] Based on the intrinsic mathematical properties of nonlinear mapping functions, this invention defines the following practical application range:

[0143] Interval 1, the micro-dynamic dominant interval, with the first sub-weight ≤ 0.05: its boundary is determined by a statistical confidence interval where the total effective sample size is greater than or equal to a preset detection threshold. Within this interval, the corresponding operation is to lock the first sub-weight assigned to the three-dimensional radiative attenuation prior tensor to a preset minimum value. High-frequency adaptive dimming commands, calculated entirely from real-time dynamic photometric trajectories, are issued via the industrial control bus to achieve the highest precision dynamic tracking of changes in local tunnel ambient brightness.

[0144] Interval two, the smooth degradation transition interval, 0.05 < first sub-weight < 0.95: its boundary is determined by the sensitive region of the exponential mapping of the scaling bias on the negative half-axis. Within this interval, the absolute value of the first derivative of the exponential function is in the peak region, characterizing a highly sensitive response to probe loss and residual expansion. The corresponding operation is to extract the continuously changing first and second sub-weights in real time and perform a weighted fusion operation. The actuator receives the mixing and smoothing instruction, relies on the automatically assigned fusion weight parameter to absorb the impact of sudden environmental changes, and eliminates the mechanical oscillations during the dynamic control of the residual distribution state and weight degradation transition.

[0145] Interval 3, the pure prior graph prediction interval, with the first sub-weight ≥ 0.95: its boundary is determined by the absolute value of the scaling deviation crossing the curvature abrupt change point, causing the intermediate exponent in the denominator of the function to approach zero, thus prompting the calculation of the first sub-weight to enter the mathematical saturation region. This mathematical saturation state is physically equivalent to the complete loss of information entropy of the probe network. The corresponding operation is to physically shield all input signal interference from the set of differential quotients of dynamic photometric trajectories, forcing the three-dimensional radiation attenuation prior tensor to take full control of the generation of the positive detection indication signal. Based on the macroscopic logical knowledge graph output, a constant safe lighting baseline is established, establishing the operational baseline under extreme blinding conditions, until the total number of effective samples again meets the statistical confidence requirements.

[0146] Figure 4This is a schematic diagram illustrating the results of a numerical verification experiment conducted in an embodiment of the present invention, targeting a nonlinear smoothing degradation weight allocation mechanism based on the spatiotemporal residual distribution matrix and the total effective sample size. This numerical verification aims to test the theoretical response characteristics of the present invention under preset standardized boundary conditions (detection threshold critical point and test vectors with superimposed multi-level meteorological distortions). Figure 4 In the diagram, the horizontal axis represents the total number of valid samples calculated, and the vertical axis represents the calculated prior weight allocation values. The dashed right-angled line and the solid smooth curve represent the data series of the first weight in the traditional Boolean scheme and the first sub-weight in this invention, respectively. Figure 4 As shown, the numerical calculation results reveal the fundamental morphological differences in the output trajectories of control commands under two different architectures. The nonlinear inverse proportional activation logic proposed in this invention can eliminate command jumps caused by rigid judgments. Specifically, within the test interval where the preset detection threshold is exceeded (in this embodiment, the horizontal coordinate approaches the critical region of 490), the calculated data series of this invention exhibits the smooth transition characteristics unique to mathematical functions (smoothly evolving from 0.050 to 0.768); in contrast, within the same region, the traditional Boolean scheme data series exhibits an absolute binary right-angle jump (instantly jumping from 0.050 to 1.000). By introducing a composite exponential mapping mechanism between the global residual norm and the negative deviation, the intrinsic properties of the mathematical function can actively absorb the impact of sudden environmental changes, reducing the step impact amplitude of control commands.

[0147] Figure 1 Using a high-information-density vector-structured hierarchical diagram, the complete computational execution path is objectively displayed, starting from the hardware perception source on the left and sequentially passing through the four core processing hubs on the right. Each main logic rectangle in the diagram represents a processing stage, and the nested sub-logic modules within it correspond to further feature-limited processes with physical constraints and mathematical mapping relationships under each main step. For example... Figure 1As shown, the starting point of the technical roadmap is the physical observation front-end at the tunnel entrance on the left. Data flows from here into the first main logic box of the roadmap, "Macroscopic Prior Tensor Inference" (corresponding to step S1). Its internal sub-modules demonstrate the acquisition of a macroscopic logical knowledge graph containing ephemeris status and meteorological cloud cover, and the calculation of the cross-scale attention weight matrix between global meteorology and local topography. This is followed by a Maxwell-constrained node feature aggregation operation on the initial three-dimensional radiation field. The second main logic box, "Microscopic Desensitization Trajectory Extraction" (corresponding to step S2), demonstrates the purification process of the microscopic data stream. This involves performing irreversible edge anonymization desensitization on the original sequence, followed by sequentially performing specular rejection verification based on pixel variation coefficients, and then normalizing the dimensions by introducing a transient velocity scalar to output an absolutely clean set of differential quotients. The third main logic box, "Bayesian Residual Assimilation" (corresponding to step S3), demonstrates... The mapping calculation of heterogeneous data is represented by a sub-module characterization system that constructs a four-dimensional spatiotemporal interpolation grid based on the intrinsic parameter matrix and constructs the prior probability density function domain by calculating the attenuation integral along the path of the virtual ray beam. Finally, the spatiotemporal residual distribution matrix of the mapping difference is obtained. The fourth main logic block, "bidirectional collaborative evolution closed loop" (corresponding to steps S4 and S5), integrates the output and iteration mechanism of this embodiment. It clearly divides two paths: the forward path compares the total number of effective samples with the preset detection threshold and uses nonlinear mapping to allocate the first sub-weight to perform smooth degradation fusion and output indication signal; the reverse path uses a high-pass filter to extract high-frequency variation components, and generates a regularization penalty term after completing the spatial clustering index verification, forcing the graph neural network to update the local aerosol scattering coefficient in the next cycle.

[0148] The computational logic involved in this application can be constructed using algorithms such as regression analysis in machine learning, establishing a mathematical model by analyzing the inherent trends and interrelationships of the collected parameters. This process can be implemented using professional computational tools (such as Python's Scikit-learn library or the R language environment). To ensure the effectiveness and accuracy of the model, techniques such as cross-validation will be used to evaluate model performance, and iterative optimization will be performed based on continuous feedback to ensure that the computational content truly reflects the inherent laws of objective data. In all computational processes, to eliminate the influence of different physical dimensions and ensure that data is compared and analyzed on the same scale, the input parameters in each formula are dimensionless. The dimensionless techniques used include, but are not limited to, min-max-normalization or Z-score standardization.

[0149] It should be emphasized that the foregoing embodiments are merely illustrative of preferred implementations of the present invention and are not intended to limit the scope of protection of the present invention. Those skilled in the art can make various modifications, equivalent substitutions, or improvements based on the technical solutions, spirit, and principles disclosed in this invention. Any such variations that do not depart from the scope of the technical solutions of this invention should fall within the protection scope of the claims of this invention. Furthermore, the technical features of the various embodiments described in the specification can be arbitrarily combined without creating technical contradictions. All possible, non-contradictory combinations should be considered within the scope of this specification; this application also provides a computer-readable storage medium storing computer program instructions thereon.

Claims

1. A method for monitoring road condition images of the black hole effect at tunnel entrances, characterized in that, The specific steps include: S1: Obtain the three-dimensional radiation attenuation prior tensor that represents the macroscopic radiation attenuation physical law within the tunnel entrance transition space defined by the target geographic coordinates and local terrain structure. S2: Obtain the original image sequence characterizing the spatiotemporal changes of photometric surface of dynamic target within the specific spatial environment; perform irreversible edge anonymization and desensitization on the original image sequence to remove identity features; and extract the set of dynamic photometric trajectory differential quotients based on the desensitized image sequence. S3: Discretize the three-dimensional radiation attenuation prior tensor into the spatiotemporal topological space corresponding to the set of differential quotients of dynamic photometric trajectories, and perform Bayesian residual assimilation calculation to generate a spatiotemporal residual distribution matrix that maps the difference between the three-dimensional radiation attenuation prior tensor and the set of differential quotients of dynamic photometric trajectories. S4: Based on the spatiotemporal residual distribution matrix, allocate the fusion weight parameters of the three-dimensional radiation attenuation prior tensor and the set of differential quotients of the dynamic photometric trajectory, and perform a weighted fusion operation based on the fusion weight parameters to generate a positive detection indication signal characterizing the global adaptive calibration absolute brightness attenuation distribution. S5: Extract the high-frequency variation components in the spatiotemporal residual distribution matrix, and use the high-frequency variation components as regularization penalty terms to generate a reverse evolution indication signal for triggering the iterative update of the local climate parameters of the three-dimensional radiation attenuation prior tensor in the next calculation cycle.

2. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 1, characterized in that: Steps to obtain the three-dimensional radiation attenuation prior tensor: Obtain macroscopic logical knowledge graph data containing target geographic coordinates, ephemeris status, and real-time meteorological cloud cover status; Using the macroscopic logical knowledge graph data as input features, a graph neural network topological reasoning process based on simplified boundary condition constraints according to Maxwell's radiative transfer equation is executed to calculate the three-dimensional radiative attenuation prior tensor with spatial depth attributes.

3. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 2, characterized in that: The steps for executing the graph neural network topological reasoning process based on simplified boundary condition constraints according to Maxwell's radiative transfer equation are as follows: extract the global meteorological feature map and the local terrain structure map from the macroscopic logical knowledge graph data; calculate the cross-scale attention weight matrix between the global meteorological feature map and the local terrain structure map; Based on the cross-scale attention weight matrix, node features are aggregated on the initial three-dimensional radiation field to generate the three-dimensional radiation attenuation prior tensor.

4. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 3, characterized in that: The steps for extracting the dynamic photometric trajectory differential quotient set based on the desensitized image sequence are: extracting the dynamic target bounding box sequence from the desensitized continuous temporal images; For each frame target region in the dynamic target bounding box sequence, calculate its monochrome variation coefficient in multiple color channels in the first color space; An arithmetic mean is performed on the monochrome coefficients of variation across multiple color channels to generate the pixel coefficients of variation for the dynamic target. Compare the pixel variation coefficient with a preset specular repulsion threshold; When the pixel variation coefficient is lower than the preset specular repulsion threshold, the apparent gray-scale temporal difference of the same dynamic target on the spatial displacement vector near the tunnel entrance is calculated to generate the set of dynamic photometric trajectory differential quotients.

5. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 4, characterized in that: The steps for calculating the apparent gray-scale temporal difference of the same dynamic target on the spatial displacement vector near the tunnel entrance are as follows: extract the motion optical flow vector of the same dynamic target between consecutive frames; calculate the transient velocity scalar based on the motion optical flow vector; introduce a dynamic normalization parameter, and perform a multiplication operation between the transient velocity scalar and the dynamic normalization parameter to generate a dimensionless velocity penalty factor. The ratio of the apparent gray-scale temporal difference to the velocity penalty factor is calculated, and the ratio is mapped to a preset normalization interval to generate the set of differential quotients of the dynamic photometric trajectory.

6. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 5, characterized in that: The steps for performing Bayesian residual assimilation calculations to generate the spatiotemporal residual distribution matrix of the mapping differences are as follows: Construct a four-dimensional spatiotemporal interpolation grid based on the preset camera intrinsic parameter transformation matrix; The three-dimensional radiation attenuation prior tensor is mapped to the prior probability density function domain of the four-dimensional spatiotemporal interpolation grid; The set of differential quotients of the dynamic photometric trajectory is converted into a likelihood observation vector; the joint posterior probability expectation deviation between the prior probability density function domain and the likelihood observation vector is solved to generate the spatiotemporal residual distribution matrix characterizing the degree of deviation between the two.

7. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 6, characterized in that: The steps for mapping the three-dimensional radiation attenuation prior tensor to the prior probability density function domain of the four-dimensional spatiotemporal interpolation grid are as follows: A virtual ray beam set is constructed with the camera's optical center as the origin; the in-path attenuation integral path of the virtual ray beam set as it penetrates the three-dimensional radiation attenuation prior tensor is calculated; based on the in-path attenuation integral path, expected values ​​and variance parameters are assigned to each time slice of the four-dimensional spatiotemporal interpolation grid to construct the prior probability density function domain.

8. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 7, characterized in that: The steps for assigning fusion weight parameters and performing weighted fusion operation are as follows: Calculate the total number of valid samples contained in the set of differential quotients of the dynamic photometric trajectory; compare the total number of valid samples with the preset detection threshold. When the total number of valid samples is greater than or equal to the preset detection threshold, a micro-dominant calculation strategy is executed: specifically, the state normalization parameters are extracted, mapped to a minimum weight value through a linear inverse proportional function and assigned to the first sub-weight corresponding to the three-dimensional radiation attenuation prior tensor, and the difference between the value 1 and the first sub-weight is assigned to the second sub-weight corresponding to the set of differential quotients of the dynamic photometric trajectory. When the total number of valid samples is lower than the preset detection threshold, a downgrade calculation strategy is executed: the first sub-weight corresponding to the three-dimensional radiation attenuation prior tensor is monotonically increased, and the second sub-weight corresponding to the set of differential quotients of dynamic photometric trajectories is monotonically decreased simultaneously to generate a positive detection indication signal.

9. The method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 8, characterized in that: The steps are as follows: Monotonically increasing the first sub-weight corresponding to the three-dimensional radiative attenuation prior tensor, and simultaneously monotonically decreasing the second sub-weight corresponding to the set of differential quotients of dynamic photometric trajectories. Obtain the total number of effective samples and the preset detection threshold; subtract the preset detection threshold from the total number of effective samples to obtain a negative deviation value with a quantitative dimension; extract the global residual norm of the spatiotemporal residual distribution matrix, and calculate and generate state normalization parameters based on the global residual norm; The product of the negative deviation value and the state normalization parameter is calculated to obtain the dimensionless scaling deviation. The scaling bias is input into a nonlinear relationship mapping function based on data fitting; The specific calculation logic for the nonlinear relational mapping function is as follows: using the natural constant as the base and the dimensionless scaling deviation as the exponent, perform an exponentiation operation to obtain an intermediate exponent result; add the value 1 to the intermediate exponent result to obtain the denominator; divide the value 1 by the denominator and extract the output quotient as the first sub-weight; perform a subtraction operation to subtract the first sub-weight from the value 1 to obtain the second sub-weight.

10. A method for monitoring road condition images of the black hole effect at tunnel entrances according to claim 9, characterized in that: The steps to generate a reverse evolution indicator signal are as follows: A time-domain high-pass filter is applied to decompose the spatiotemporal residual distribution matrix to isolate the high-frequency variation components with frequencies higher than the first cutoff frequency; Calculate the spatial clustering index of the high-frequency variation components; When the spatial clustering index meets the preset convergence condition, the average absolute value of the residual of the high-frequency variation component is obtained; a loss mapping transformation constant is introduced to eliminate the difference in physical dimensions, and the spatial clustering index, the average absolute value of the residual, and the loss mapping transformation constant are multiplied together to generate a dimensionless pure numerical result. The dimensionless pure numerical result is established as the regularization penalty term, and the regularization penalty term is injected into the physical loss function of the graph neural network in an additive form. Then, the weight gradient of the hidden layer nodes of the network is calculated through the error backpropagation algorithm, thereby forcing the graph neural network to update the local aerosol scattering coefficient, which is represented as the network node parameter, in the next calculation cycle. The steps for calculating the spatial clustering index of the high-frequency variation components are as follows: Determine the three-dimensional Euclidean distance matrix of the high-frequency mutation component in the three-dimensional spatial coordinate system; for each target node in the high-frequency mutation component, based on the three-dimensional Euclidean distance matrix, count the number of adjacent nodes whose distance is less than the preset search radius; When the number of adjacent nodes is greater than a preset core density threshold, the target node is marked as a core node; Summarize the total number of the marked core nodes; By introducing a positive real number safety cushion parameter, the total number of nodes of the high-frequency variation component is added to the positive real number safety cushion parameter to obtain a safety denominator that is forced to be non-zero. The ratio of the total number of core nodes to the safety denominator is calculated and used as the spatial clustering index.

Images