Method for real-time monitoring of material blending uniformity based on laser scattering

Abstract

The application discloses a material blending uniformity real-time monitoring method based on laser scattering and belongs to the technical field of material blending uniformity real-time monitoring. The method comprises the following steps: collecting a multi-angle laser scattering signal set after a laser beam penetrates through a material to be detected; pre-processing the multi-angle laser scattering signal set to generate a multispectral fusion feature map; inputting the multispectral fusion feature map into a chaotic attractor reconstructor to extract an orbit distribution topological structure of a material particle group; performing singular spectrum analysis on the orbit distribution topological structure to locate a local component concentration abnormal area in a material blending process; calculating statistical self-similarity parameters of the local component concentration abnormal area to construct a dimensionless index representing overall material blending uniformity. The method can reveal chaotic mixing dynamics characteristics of the particle group and realize accurate evaluation of the uniformity based on abnormal area self-similarity.

Description

Technical Field

[0001] This invention relates to the field of real-time monitoring technology for material blending uniformity, specifically a method for real-time monitoring of material blending uniformity based on laser scattering. Background Technology

[0002] The uniformity of material blending directly affects the performance consistency of multi-component particulate products, making its real-time monitoring crucial in chemical, pharmaceutical, and food processing industries. Existing online monitoring technologies generally employ single-angle or finite discrete-angle laser scattering methods. These methods determine the mixing endpoint by detecting the time-series fluctuations in scattered light intensity from a fixed direction and calculating statistical indicators such as the coefficient of variation. Single-angle detection schemes assume a simple linear relationship between scattered light intensity and component concentration. However, in real-world systems with wide particle size distributions and irregular shapes, the selectivity of a single-angle light intensity for local particle aggregation structures is limited, resulting in the loss of a significant amount of spatial structural information contained within the scattering angle distribution. Multi-angle detection techniques typically analyze each signal separately and then average or weighted fuse them, failing to uncover the coupling relationships between signals from different angles and the implicit nonlinear dynamics of particle swarms. The convection and diffusion processes of particles within mixing equipment are essentially chaotic dynamic systems. Local component concentration anomalies often appear as isolated rich clusters or stripes, exhibiting statistical self-similarity in their spatial distribution. Current technologies lack effective means to reconstruct the phase space trajectory of this system from multi-angle scattering sequences, thus making it impossible to physically locate these anomalous regions and extract their geometric structural features. Conventional uniformity indices, such as the coefficient of variation, only reflect the second-order statistics of global concentration deviation. In cases of localized high aggregation where the overall mean remains close to the target value, they may output false acceptable signals, failing to accurately assess mixing quality. Existing methods cannot simultaneously reconstruct the orbital topology of chaotic particle swarm mixing from multi-angle scattering information, nor can they quantify the spatial self-similarity of local anomalous regions into an overall uniformity index, resulting in significant deficiencies in the sensitivity and accuracy of real-time monitoring. Summary of the Invention

[0003] This invention provides a real-time monitoring method for material mixing uniformity based on laser scattering, aiming to solve the problems of existing technologies, such as difficulty in reconstructing the nonlinear dynamic topological characteristics of particle group mixing process from multi-angle laser scattering signals, difficulty in accurately locating local component concentration abnormal regions, and lack of overall uniformity index that can comprehensively reflect the self-similarity of abnormal space.

[0004] To achieve the above objectives, this invention provides the following technical solution: This invention proposes a real-time monitoring method for material blending uniformity based on laser scattering. It collects a set of multi-angle laser scattering signals after a laser beam penetrates the material under test, and preprocesses this set to generate a multispectral fusion feature map. The multispectral fusion feature map is input into a chaotic attractor reconstructor to extract the orbital distribution topology of the material particle group. Singular spectrum analysis is performed on the orbital distribution topology to locate local component concentration anomalies during the material blending process. Statistical self-similarity parameters for these local component concentration anomalies are calculated to construct a dimensionless exponent representing the overall blending uniformity of the material. Through this method, online, continuous, and interference-free monitoring of material blending uniformity can be achieved, significantly improving the ability to distinguish microscopic differences in component distribution during blending, and overcoming the shortcomings of traditional sampling and detection methods, such as lag, sample destruction, and inability to reflect the overall picture.

[0005] As a preferred embodiment of the present invention, the process of preprocessing a multi-angle laser scattering signal set to generate a multispectral fusion feature map includes: segmenting the laser scattering signal at each angle into a time series to obtain scattering intensity subsequences with equal time intervals; extracting statistical features from the scattering intensity subsequences, including peak intensity, full width at half maximum (FWHM), and skewness coefficient; encoding the statistical features into a two-dimensional matrix according to spatial angular order to obtain an initial feature map; and performing three-channel pseudo-color mapping on the initial feature map, where the first channel corresponds to peak intensity, the second channel corresponds to FWHM, and the third channel corresponds to skewness coefficient, thereby generating a multispectral fusion feature map. This processing method compresses multidimensional scattering information into an image representation that integrates spatial and intensity distribution characteristics, greatly enhancing the differential feature representation of scattering fingerprints of different component particle groups, thus enabling higher sensitivity and stability in subsequent dynamic evolution analysis of the mixing state.

[0006] As a preferred embodiment of the present invention, the process of extracting the orbital distribution topology of a particle swarm from a chaotic attractor by inputting a multispectral fusion feature map into the chaotic attractor reconstructor includes pixel-level normalization of the multispectral fusion feature map to generate a normalized feature map. The normalized feature maps are arranged in chronological order to form a high-dimensional phase space vector sequence. A delayed embedding transformation is performed on the high-dimensional phase space vector sequence to reconstruct the phase space trajectory of the nonlinear dynamic system. The neighborhood connection matrix of the phase space trajectory is extracted, which records the connection relationship between each phase space point and its nearest neighbor. The neighborhood connection matrix is ​​projected onto a low-dimensional manifold to obtain the orbital distribution topology. This method extracts the structural features of the chaotic attractor from the disordered scattering signal, mapping the discrete motion of the particle swarm to a continuous topological network. It can reveal the inherent dynamic organization of the mixed system in a geometrically intuitive way and sensitively capture the microscopic inhomogeneities hidden beneath the macroscopic homogeneous appearance.

[0007] As a preferred technical solution of the present invention, the process of locating local component concentration anomaly regions by performing singular spectral analysis on the orbital distribution topology includes unfolding the orbital distribution topology into a one-dimensional time series, which contains topological feature values ​​arranged in the order of acquisition. A trajectory matrix is ​​constructed from this one-dimensional time series, with each row element composed of subsequences truncated by a sliding window. Singular value decomposition is performed on the trajectory matrix to generate a singular value sequence and corresponding left and right singular vectors. Trend components and detail components are reconstructed based on the left and right singular vectors; the trend component reflects the uniformity baseline, and the detail component reflects local fluctuations. Spatiotemporal coordinate points where the amplitude of the detail component exceeds a threshold are marked as local component concentration anomaly regions. This step can automatically extract the baseline representing the main mixing trend and the fluctuation component representing local anomalies from the global signal, achieving accurate spatiotemporal location of concentration defect regions such as minor segregation and aggregation, significantly improving the automation level and reliability of anomaly identification.

[0008] As a preferred technical solution of the present invention, the process of calculating the statistical self-similarity parameters of local component concentration anomaly regions and constructing a dimensionless index includes binarizing each local component concentration anomaly region to generate a binary anomaly region mask. The binary anomaly region mask is covered by grids of different scales, and the minimum number of grids covering the binary anomaly region at each grid scale is recorded. A linear regression is performed with the logarithm of the grid scale as the independent variable and the logarithm of the minimum number of grids as the dependent variable, and the slope of the fitted line is calculated as the fractal dimension. This fractal dimension is then compared with the fractal dimension of a standard homogeneous material to generate a dimensionless index representing the overall blending homogeneity of the material. Based on the geometric irregularity and self-similarity characteristics of the anomaly region, this index can quantify the global blending homogeneity state in a comprehensive numerical value, providing a unified, objective, and traceable evaluation standard for the blending quality of different batches and at different times.

[0009] As a preferred embodiment of the present invention, the method for acquiring a multi-angle laser scattering signal set includes setting up a ring-array laser receiver, which comprises multiple photodetectors distributed at equal angular intervals. A monochromatic laser source is excited and incident on the material under test along a vertical direction, causing the material to scatter light. The ring-array laser receiver is synchronously triggered, causing each photodetector to simultaneously acquire scattering intensity data at the corresponding angle. The scattering intensity data is encapsulated according to the angle index and acquisition timestamp to generate a multi-angle laser scattering signal set. This acquisition architecture achieves multi-directional synchronous capture of the material's scattered light field, ensuring consistency on the time axis and integrity of spatial information, providing a solid data foundation for subsequent high-fidelity spatial data construction.

[0010] As a preferred embodiment of the present invention, before inputting the multispectral fusion feature map into the chaotic attractor reconstructor, the method further includes blind source separation processing of the multispectral fusion feature map to remove common-mode noise components introduced by ambient stray light. Median filtering is performed on abnormal pixels in the multispectral fusion feature map to repair signal loss areas caused by dust or air bubbles. Spatial resolution downsampling is performed on the multispectral fusion feature map to obtain a uniform-sized input feature map. These preprocessing operations effectively suppress artifacts caused by environmental interference and equipment defects, improve the purity and consistency of the input data, and ensure the accuracy and robustness of the chaotic attractor reconstruction.

[0011] As a preferred embodiment of the present invention, after extracting the orbital distribution topology of the material particle group, the method further includes calculating the global Lyapunov exponent of the orbital distribution topology, which characterizes the degree of chaos in the material mixing system. Isolated nodes in the neighborhood connectivity matrix of the orbital distribution topology are removed to obtain a simplified topology. Spectral clustering analysis is performed on the simplified topology to divide the material particle group into multiple sub-clusters with similar dynamic behavioral characteristics. This additional analysis not only allows for macroscopic determination of whether the system is in an ordered or chaotic mixing state, but also identifies particle groups with similar motion patterns at the microscopic level, providing a direct basis for diagnosing mixing dead zones or components with poor following behavior.

[0012] As a further optimization of the above technical solution, the process of performing spectral clustering analysis on the simplified topology includes constructing a similarity matrix of the simplified topology, where each element of the similarity matrix is ​​an exponential decay function value of the topological distance between any two nodes. The Laplace matrix of the similarity matrix is ​​calculated, and eigenvalue decomposition is performed on this Laplace matrix to extract the eigenvectors corresponding to the k smallest eigenvalues. These k eigenvectors are arranged row-wise to form an eigenvector matrix, and k-means clustering is performed on each row of the eigenvector matrix to obtain the partitioning results of multiple sub-clusters. This algorithm can dynamically segment particle groups while preserving topological relationships, and the partitioned sub-clusters accurately reflect the different hybrid dynamic regions within the material.

[0013] As a further preferred embodiment of the present invention, the aforementioned blind source separation processing employs an independent component analysis algorithm, and the number of independent components is set to three. In this way, the overlapping useful signals and environmental noise in the originally acquired multi-channel scattering signals can be effectively separated, further improving the quality of the feature map.

[0014] The technical effects and advantages provided by the present invention in the above technical solution are as follows: The multi-angle laser scattering signal set is segmented according to time series, and statistical features such as peak intensity, full width at half maximum (FWHM), and skewness coefficient are extracted and encoded in spatial angular order to generate a multispectral fusion feature map. Pixel-level normalization is performed on this feature map, and a high-dimensional phase space vector sequence is formed in time series. The phase space trajectory of the nonlinear dynamic system is reconstructed through delayed embedding transformation. The neighborhood connectivity matrix formed by each phase point and its nearest neighbor is extracted and projected onto a low-dimensional manifold to obtain the orbital distribution topology of the particle swarm. This processing explicitly maps the geometric features of the chaotic attractors hidden in the scattering signal into an analyzable topological network. Subtle changes in the collective dynamic state of the particle swarm can cause local abrupt changes in the topological connectivity or the generation of isolated nodes, enabling sensitive indication of local component concentration anomalies with high spatial resolution. Compared to traditional methods that rely on statistical thresholds, this topological structure is insensitive to illumination fluctuations and material color changes, significantly enhancing the robustness of anomaly detection. The orbital distribution topology is unfolded into a one-dimensional time series, and a trajectory matrix is ​​constructed. Singular value decomposition is performed to obtain trend components reflecting the uniformity baseline and detail components characterizing local fluctuations. Local anomaly regions in component concentration are precisely marked based on whether the amplitude of the detail components exceeds a threshold. Singular spectral analysis adaptively removes system baseline drift, ensuring anomaly localization is unaffected by gradual changes in overall concentration. Based on multi-scale grid coverage of a binary mask over the anomaly region, the rate of change of the logarithm of the required grid number with the logarithm of the scale is calculated as the fractal dimension. This is then compared with the fractal dimension of the standard homogeneous material to generate a dimensionless uniformity index. This index is rooted in the spatial self-similarity of local anomaly regions and is highly sensitive to inhomogeneous structures such as clustering and layered separation, effectively detecting hidden defects that traditional coefficients of variation cannot identify. Simultaneously, its dimensionless nature eliminates the influence of laser power, material optical properties, and differences in mixing equipment, making the uniformity obtained under different process conditions directly comparable. Attached Figure Description

[0015] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments recorded in this invention. For those skilled in the art, other drawings can be obtained based on these drawings.

[0016] Figure 1 This is a flowchart of a method for real-time monitoring of material mixing uniformity based on laser scattering; Figure 2 This is a schematic diagram of a ring array laser scattering detection system; Figure 3 This is a flowchart for generating multispectral fusion feature maps based on multi-angle laser scattering signals; Figure 4 This is a flowchart of the multispectral fusion characteristic phase space reconstruction and particle group sub-cluster division process in the material blending process; Figure 5 This is a flowchart of a method for detecting abnormal concentrations of local components during the material blending process. Figure 6 This is a flowchart for calculating the dimensionless index of material blending uniformity. Detailed Implementation

[0017] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0018] See Figure 1 This invention provides a real-time monitoring method for material blending uniformity based on laser scattering, comprising: acquiring a set of multi-angle laser scattering signals after a laser beam penetrates the material to be tested; preprocessing the set of multi-angle laser scattering signals to generate a multispectral fusion feature map; inputting the multispectral fusion feature map into a chaotic attractor reconstructor to extract the orbital distribution topology of the material particle group; performing singular spectrum analysis on the orbital distribution topology to locate local component concentration anomaly regions during the material blending process; calculating the statistical self-similarity parameters of the local component concentration anomaly regions to construct a dimensionless index representing the overall blending uniformity of the material.

[0019] Example 1: In specific implementation, please refer to Figure 2 A ring-shaped laser receiver is configured. This receiver comprises N photodetectors, distributed at equal angular intervals on the same circumference, where N is an integer greater than 1, and the angular interval is (360 / N) degrees. Each photodetector has a unique angle index, numbered sequentially from 0 to N-1. The photodetector with index i corresponds to a detection azimuth angle of i × (360 / N) degrees, where i is an integer from 0 to N-1. The circumferential plane of the ring-shaped laser receiver is perpendicular to the conveying direction of the material being measured, and the center of the ring-shaped laser receiver is located on the incident optical axis of the monochromatic laser source.

[0020] In practice, a monochromatic laser source is excited and incident on the material to be tested in a vertical direction. The monochromatic laser source outputs a continuous laser beam with a fixed wavelength, and the beam diameter is smaller than the cross-sectional dimension of the material to be tested in the direction perpendicular to the beam. The distance between the exit port of the monochromatic laser source and the upper surface of the material to be tested is determined according to the Rayleigh length of the laser beam and the working distance of the material to be tested, so that the laser beam is at the beam waist position when it is incident on the material to be tested. The material to be tested moves at a constant speed on the conveying device. After the monochromatic laser source is excited, the laser beam penetrates the material to be tested. Different component particles inside the material to be tested scatter the laser beam, and the scattered light propagates outward from the material to be tested at different angles.

[0021] In practical implementation, a synchronous triggering ring array laser receiver is used. A synchronous triggering controller is set up, which simultaneously outputs trigger pulse signals to both the monochromatic laser source and the ring array laser receiver. At the rising or falling edge of each trigger pulse, the monochromatic laser source outputs a laser pulse, and the N photodetectors in the ring array laser receiver synchronously initiate a scattering intensity acquisition. Each photodetector integrates a photoelectric conversion element, an amplifier, and an analog-to-digital converter. The photoelectric conversion element converts the scattered light signal at the corresponding angle into a current signal, the amplifier amplifies the current signal and converts it into a voltage signal, and the analog-to-digital converter converts the voltage signal into digital scattering intensity data. The repetition frequency of the trigger pulse output by the synchronous triggering controller is set to f Hz based on the speed of the material being measured and the required monitoring time resolution, and the duty cycle of the trigger pulse is 50%.

[0022] In practice, the scattering intensity data is encapsulated according to angle index and acquisition timestamp to generate a multi-angle laser scattering signal set. After each synchronous trigger acquisition, the ring array laser receiver outputs N scattering intensity data, each with a corresponding angle index. Simultaneously, the synchronous trigger controller records the timestamp of this acquisition. The N scattering intensity data, their corresponding angle indices, and the acquisition timestamp are combined to form an acquisition record. Multiple acquisition records formed from multiple acquisitions are arranged in ascending order of timestamps, generating a multi-angle laser scattering signal set with timestamps as rows and angle indices as columns. The multi-angle laser scattering signal set is stored in matrix form, with the number of rows being the total number of acquisitions T and the number of columns being N. The element S(t,b) in the matrix represents the scattering intensity data of the t-th acquisition with angle index b, where t ranges from 1 to T, and b ranges from 0 to N-1.

[0023] Example 2: In specific implementation, please refer to Figure 3The laser scattering signal at each angle is segmented into time series to obtain scattering intensity subsequences with equal time intervals. The scattering intensity data column with angle index b is extracted from the multi-angle laser scattering signal set, where b ranges from 0 to N-1. Let the total number of acquisitions be T, and the scattering intensity data column with angle index b contains T scattering intensity data points arranged in order of acquisition timestamp. A subsequence length parameter L is set, where L is an integer greater than 1 and less than T. The subsequence length parameter L is determined based on the ratio of the velocity of the material to be measured to the diameter of the laser beam, ensuring that the material sampling volume corresponding to each scattering intensity subsequence is approximately equal. The scattering intensity data column with angle index b is segmented without overlap using a step size of L, obtaining M scattering intensity subsequences with equal time intervals, where M is equal to T divided by L and rounded down. The p-th scattering intensity subsequence contains L scattering intensity data points, where p ranges from 1 to M. The above time series segmentation operation is performed on N angle indices respectively, resulting in a total of N multiplied by M scattering intensity subsequences.

[0024] In practice, statistical features of the scattering intensity subsequences are extracted, including peak intensity, full width at half maximum (FWHM), and skewness coefficient. The maximum value among the L scattering intensity data points in the p-th scattering intensity subsequence is searched and taken as the peak intensity of the p-th scattering intensity subsequence. Two time positions corresponding to 50% of the peak intensity value are searched in the p-th scattering intensity subsequence; these two time positions are the left and right half-width positions, respectively. The difference between the time coordinates corresponding to the right and left half-width positions is taken as the FWHM of the p-th scattering intensity subsequence. The skewness coefficient is calculated for the p-th scattering intensity subsequence, expressed by the following formula:

[0025]

[0026] in, This represents the skewness coefficient of the p-th scattering intensity subsequence at angle index b. This represents the subsequence length parameter. This represents the q-th scattering intensity data in the p-th scattering intensity subsequence at angle index b, where q ranges from 1 to L. This represents the arithmetic mean of the L scattering intensity data in the p-th scattering intensity subsequence at angle index b. For each angle index b and each scattering intensity subsequence p, peak intensity extraction, full width at half maximum (FWHM) extraction, and skewness coefficient calculation are performed to obtain N multiplied by M sets of statistical features. Each set of statistical features contains a peak intensity value, a FWHM value, and a skewness coefficient value.

[0027] In the specific implementation, the statistical features are encoded into a two-dimensional matrix according to spatial angular order to obtain the initial feature map. A two-dimensional matrix with N rows and 3M columns is constructed, where the b-th row of the two-dimensional matrix corresponds to the scattering direction with angle index b. The peak intensity value of the scattering intensity subsequence with angle index b and subsequence number p is assigned to the element in the b-th row and (3p-2)-th column of the two-dimensional matrix; the full width at half maximum (FWHM) value of the scattering intensity subsequence with angle index b and subsequence number p is assigned to the element in the b-th row and (3p-1)-th column of the two-dimensional matrix; and the skewness coefficient value of the scattering intensity subsequence with angle index b and subsequence number p is assigned to the element in the b-th row and 3p-th column of the two-dimensional matrix. The value of b ranges from 0 to N-1, and the value of p ranges from 1 to M. The row order of the two-dimensional matrix is ​​consistent with the increasing order of the spatial angles, and the column order of the two-dimensional matrix is ​​arranged according to the time order of the subsequences. Each group contains three columns of eigenvalues: peak intensity, full width at half maximum (FWHM), and skewness coefficient. This two-dimensional matrix is ​​the initial feature map.

[0028] In the specific implementation, a three-channel pseudo-color mapping is performed on the initial feature map, where the first channel corresponds to peak intensity, the second channel corresponds to full width at half maximum (FWHM), and the third channel corresponds to skewness coefficient, generating a multispectral fusion feature map. Three feature planes are separated from the initial feature map, each composed of peak intensity feature values, FWHM feature values, and skewness coefficient feature values. The peak intensity feature plane consists of all (3p-2)th columns of the initial feature map, where p ranges from 1 to M. The FWHM feature plane consists of all (3p-1)th columns of the initial feature map, and the skewness coefficient feature plane consists of all (3p)th columns of the initial feature map. Linear normalization is performed on the peak intensity feature plane, mapping each element value to an integer range of 0 to 255, resulting in the first channel grayscale image. Similarly, linear normalization is performed on the FWHM feature plane, mapping each element value to an integer range of 0 to 255, resulting in the second channel grayscale image. The skewness coefficient feature plane is linearly normalized, mapping each element value in the skewness coefficient feature plane to an integer range of 0 to 255, resulting in the third channel grayscale image. The first, second, and third channel grayscale images are then superimposed in channel order to generate a three-channel pseudo-color image, which is the multispectral fusion feature map.

[0029] In the specific implementation, before inputting the multispectral fusion feature map into the chaotic attractor reconstructor, blind source separation processing is performed on the multispectral fusion feature map. An independent component analysis algorithm is used, with the number of independent components set to 3, to remove common-mode noise components introduced by ambient stray light. The three channel images of the multispectral fusion feature map are each expanded row-wise into three one-dimensional observation signal vectors, each with a length of N multiplied by M. These three one-dimensional observation signal vectors are combined into a 3-row, N multiplied by M-column observation signal matrix, where each row corresponds to a signal sequence from one channel. The observation signal matrix is ​​then centered to ensure that the mean of each row is zero, resulting in a centered matrix. The centered matrix is ​​then whitened, its covariance matrix is ​​calculated, eigenvalue decomposition is performed, and the whitening matrix is ​​calculated using the eigenvalues ​​and eigenvectors. The whitening matrix is ​​then multiplied by the centered matrix to obtain the whitened matrix. An independent component analysis (ICA) algorithm based on maximizing negative entropy is employed, with the number of independent components set to three. This number is chosen because the multispectral fusion feature map has three channels, and the desired independent signal sources include two components related to material scattering characteristics and one component related to ambient stray light. The unmixing matrix is ​​solved, with three rows and three columns. Each row vector of the unmixing matrix is ​​iteratively optimized to maximize the non-Gaussianity metric of the separated independent component vectors. The unmixing matrix is ​​multiplied by the whitened matrix to obtain a 3-row, N-by-M-column ICA matrix. The power spectral density of each independent component vector in the ICA matrix is ​​calculated. Independent component vectors with significant peaks in power spectral density around 50 Hz or 60 Hz are identified as common-mode noise components introduced by ambient stray light, and all element values ​​in the independent component vectors corresponding to these common-mode noise components are set to zero. The inverse of the unmixing matrix is ​​used to remap the processed independent component matrix back to the original observation space, resulting in a denoised observation signal matrix. Each row of the denoised observation signal matrix is ​​then rearranged into an N-row, M-column two-dimensional image, restoring the denoised three-channel feature map.

[0030] In practice, median filtering is performed on abnormal pixels in the multispectral fusion feature map to repair signal loss areas caused by dust or air bubbles. Abnormal pixel detection is performed on each channel of the three-channel feature map after blind source separation. A 3x3 processing window is taken, centered on the current pixel to be detected. The median of the pixel values ​​of the remaining 8 pixels (excluding the center pixel) within the processing window is calculated. If the absolute value of the difference between the center pixel's value and the median is greater than a preset anomaly threshold, the center pixel is marked as an abnormal pixel. The preset anomaly threshold is set to three times the standard deviation of all pixel values ​​in that channel image. The pixel values ​​marked as abnormal are replaced with the median of all non-abnormal pixel values ​​within the processing window. This abnormal pixel detection and replacement operation is performed on every pixel in the three channel images until all pixels in all channel images have been processed.

[0031] In the specific implementation, the multispectral fusion feature map is spatially downsampled to obtain a uniform-sized input feature map. The height of the target feature map is set to H pixels and the width to W pixels, where H and W are both positive integers. For each channel of the three-channel feature map after median filtering and substitution, bicubic interpolation is performed to scale each channel image from its original size of N rows and M columns to H rows and W columns. The bicubic interpolation determines the reference pixel position by mapping the target pixel coordinates back to the coordinate system of the original channel image. The pixel value of the target pixel at coordinates (h, w) is obtained by weighted summation of the pixel values ​​of the 16 surrounding pixels (4 rows and 4 columns) in the original channel image at the corresponding position using a cubic polynomial. After bicubic interpolation of the three channel images, they are combined into a three-channel H-row W-column input feature map. The three channels of the input feature map correspond to the processed peak intensity correlation component, full width at half maximum (FWHM) correlation component, and skewness coefficient correlation component, respectively.

[0032] Example 3: In specific implementation, please refer to Figure 4 The multispectral fusion feature map is subjected to pixel-level normalization to generate a normalized feature map. The multispectral fusion feature map contains three channels, each with H rows and W columns of pixels. Pixel-level normalization is performed on each channel image of the multispectral fusion feature map, mapping the pixel value of each pixel in each channel image from its original numerical range to the closed interval [0,1]. Let the channel number be c, and the value of c be one of the first, second, and third channels. The original pixel value of the pixel located in the u-th row and v-th column of the channel image with channel number c is denoted as... Where u is an integer ranging from 1 to H, and v is an integer ranging from 1 to W. Retrieve the maximum pixel value of all pixels in the image of channel c. and minimum pixel value ,in, It is the maximum value of H multiplied by W original pixel values ​​in the image of channel c. This is the minimum value of H multiplied by W original pixel values ​​in the channel image with channel number c. Normalization is performed on each pixel in the channel image with channel number c, and the normalized pixel value is... It was obtained in the following way: After pixel-level normalization is performed on the three channel images respectively, the three normalized channel images are combined into a normalized feature map with a size of H rows, W columns, and 3 channels.

[0033] In practice, the normalized feature maps are arranged in chronological order to form a high-dimensional phase space vector sequence. During the material blending process, multispectral fusion feature maps are continuously acquired for R time frames, where R is an integer greater than 1. Pixel-level normalization is performed on each of the R time frames' multispectral fusion feature maps to obtain R normalized feature maps. The R normalized feature maps are arranged in chronological order of their acquisition timestamps. The time frame number corresponding to the r-th normalized feature map is r, where r ranges from 1 to R. For each time frame, the normalized feature map is vectorized. The H rows and W columns of the r-th normalized feature map, consisting of 3 channels (H x W x 3), are expanded into a one-dimensional vector in a fixed order. This fixed order is: first, traversing the first, second, and third channels in the channel dimension; then, traversing from column 1 to column W in the column dimension; and finally, traversing from row 1 to row H in the row dimension. This yields the r-th normalized feature vector, with a length of D, where D equals H x W x 3. These R normalized feature vectors are then arranged in ascending order of time frame number r to form a high-dimensional phase space vector sequence, which contains R D-dimensional vectors.

[0034] In practice, a delayed embedding transformation is performed on the high-dimensional phase space vector sequence to reconstruct the phase space trajectory of the nonlinear dynamic system. An embedding dimension parameter d is set, where d is an integer greater than 1 and less than R, representing the number of dimensions in the reconstructed phase space. A delay time parameter τ is set, where τ is a positive integer, representing the number of vectors skipped when selecting adjacent coordinate components in the high-dimensional phase space vector sequence. The delayed embedding transformation maps the high-dimensional phase space vector sequence to a phase space trajectory matrix. The r-th row of the phase space trajectory matrix is ​​a row vector of dimension d multiplied by D. The r-th row of the phase space trajectory matrix is ​​formed by concatenating the r-th normalized eigenvector, the (r+τ)-th normalized eigenvector, the (r+2τ)-th normalized eigenvector, up to the (r+(d-1)τ)-th normalized eigenvector, where r ranges from 1 to [R-(d-1)τ], and the number of rows in the phase space trajectory matrix is ​​[R-(d-1)τ]. The delay time parameter τ is determined using the mutual information method. The mutual information value between two normalized eigenvectors with a time interval of τ in the high-dimensional phase space vector sequence is calculated, and a curve showing the mutual information value changing with τ is plotted. The value of τ corresponding to the first local minimum of the mutual information value is taken as the value of the delay time parameter τ. The embedding dimension parameter d is determined using the pseudo-nearest neighbor method. As the embedding dimension gradually increases from 2, the rate of change of the distance between the normalized eigenvector and its nearest neighbor in the reconstructed phase space is calculated for each embedding dimension. When the proportion of nearest neighbors whose distance change rate exceeds a preset threshold decreases to below a preset proportion threshold, the corresponding embedding dimension is the value of the embedding dimension parameter d. The preset distance change rate threshold is set to 10, and the preset proportion threshold is set to 5%.

[0035] In practical implementation, a neighborhood connection matrix is ​​extracted from the phase space trajectory. This matrix records the connection relationship between each phase space point and its nearest neighbor. Each row of the phase space trajectory matrix represents a phase space point, and the total number of phase space points is S, where S equals [R-(d-1)τ]. An S-row, S-column neighborhood connection matrix is ​​constructed, and the elements of the neighborhood connection matrix... Let represent the connection between the a-th phase space point and the b-th phase space point, where a and b are integers from 1 to S. For each phase space point a, calculate the Euclidean distance between phase space point a and the remaining S-1 phase space points in the phase space trajectory matrix. The Euclidean distance is calculated as the L2 norm of the difference between two row vectors. Search for the minimum value among the S-1 Euclidean distances, and denote the number of the phase space point that produces the minimum value as . , Let S be integers from 1 to S, not equal to a. In the neighborhood connectivity matrix, the a-th row and the... Column elements Set the value to 1, and set the elements of the remaining columns in row a to 0. Perform the above operation on all S phase space points to complete the construction of the neighborhood connectivity matrix.

[0036] In practice, the neighborhood connectivity matrix is ​​projected onto a low-dimensional manifold to obtain the orbital distribution topology. A local linear embedding algorithm is used to perform dimensionality reduction projection on the neighborhood connectivity matrix. The low-dimensional manifold dimension parameter k is set, where k is a positive integer greater than 1 and less than S, representing the spatial dimension after projection. For each phase space point a in the neighborhood connectivity matrix, the nearest neighbor of phase space point a is found based on the connectivity relationships in the a-th row of the neighborhood connectivity matrix. Then find the nearest neighbor point in the phase space trajectory matrix. The k-1 phase space points with the smallest Euclidean distance, and the phase space point a and its nearest neighbor. The k neighboring points of phase space point a are collectively formed. A linear reconstruction of phase space point a is performed using these k neighboring points. The linear reconstruction weight coefficients are calculated by minimizing the reconstruction error while satisfying the constraint that the sum of the weight coefficients is 1. During the calculation, a local covariance matrix is ​​constructed for each neighboring point set, and a regularization term is introduced. The value of the regularization term is set to 0.001 times the trace of the local covariance matrix. The linear reconstruction weight coefficients of the S phase space points are combined into an S-row, S-column weight matrix, where the elements corresponding to non-neighboring points are 0. A sparse symmetric matrix M is constructed based on the weight matrix. The sparse symmetric matrix M is an S-row, S-column square matrix, calculated as the identity matrix minus the transpose of the weight matrix multiplied by the identity matrix minus the weight matrix. Eigenvalue decomposition is performed on the sparse symmetric matrix M, extracting the eigenvectors corresponding to the 2nd to (k+1)th eigenvalues ​​from the smallest to the largest, resulting in k eigenvectors, each of length S. Arrange the k eigenvectors into an S-row, k-column matrix. This S-row, k-column matrix is ​​the low-dimensional manifold projection result. Each row of the low-dimensional manifold projection result is a k-dimensional vector, representing the coordinates of a phase space point in the low-dimensional manifold. The set of all S k-dimensional vectors is the orbital distribution topology.

[0037] In practice, the global Lyapunov exponent of the orbital distribution topology is calculated, which characterizes the degree of chaos in the material mixing system. Reference points are selected from the S k-dimensional vectors corresponding to the orbital distribution topology, with the number of reference points being an integer value rounded down from one-tenth of S. For each reference point, the nearest neighbor is searched within the orbital distribution topology. The time interval between the nearest neighbor and the reference point is greater than the average period, which is determined by calculating the first zero-crossing point of the autocorrelation function of the distance between adjacent time points in the orbital distribution topology. The divergence of the distance between the reference point and its nearest neighbor is tracked after several evolutionary steps, with the evolutionary step size set to an integer increment from 1 to a preset maximum evolutionary step size. For each evolutionary step size, the average logarithmic distance between all reference points and their nearest neighbors at that evolutionary step size is calculated. Using the evolution step size as the x-axis and the mean logarithmic distance as the y-axis, within the linear growth range of the mean logarithmic distance as the evolution step size changes, the slope of the fitted line is obtained by least squares linear fitting, and this slope is the estimated value of the global Lyapunov exponent.

[0038] In practice, isolated nodes are removed from the neighborhood connectivity matrix of the track distribution topology to obtain a simplified topology. Each phase space point number in the neighborhood connectivity matrix is ​​traversed, and the sum of all elements in row a and column a is checked to see if they are both zero. If both the sum of elements in row a and column a are zero, the phase space point number 'a' is marked as an isolated node. In the S x k matrix corresponding to the track distribution topology, the row vectors corresponding to the phase space point numbers marked as isolated nodes are removed. After removing all row vectors corresponding to isolated nodes, the remaining S' row vectors form an S' x k matrix, where S' equals S minus the number of isolated nodes. This S' x k matrix represents the simplified topology.

[0039] In practical implementation, spectral clustering analysis is performed on the simplified topology to divide the material particle group into multiple sub-clusters with similar dynamic behavior characteristics. A similarity matrix of the simplified topology is constructed, which is an S' x S' matrix. The elements of the similarity matrix... This represents the similarity between the e-th node and the f-th node in the simplified topology, where e and f are both integers from 1 to S'. The calculation uses an exponential decay function, expressed as:

[0040]

[0041] in, Let represent the k-dimensional coordinate vector of the e-th node in the simplified topology. Let f represent the k-dimensional coordinate vector of the f-th node in the simplified topology. express and The Euclidean distance between them Gaussian kernel bandwidth parameters Set as the median of the Euclidean distances between all pairs of nodes in the simplified topology. Diagonal elements of the similarity matrix. Set all to 0.

[0042] In practice, the Laplacian matrix of the similarity matrix is ​​calculated, and eigenvalue decomposition is performed on the Laplacian matrix to extract the eigenvectors corresponding to the first c smallest eigenvalues. c is the preset number of sub-clusters, and its value is determined based on the known number of component types in the material sample. c is an integer greater than 1 and less than S'. The degree matrix is ​​calculated; it is a diagonal matrix of S' rows and S' columns. The diagonal elements of the degree matrix... It equals the sum of the values ​​of all elements in the e-th row of the similarity matrix, that is... The degree matrix has all off-diagonal elements of 0. The non-normalized Laplacian matrix L is calculated as the degree matrix minus the similarity matrix, i.e., Perform eigenvalue decomposition on the non-normalized Laplacian matrix L to obtain S' eigenvalues ​​and their corresponding S' eigenvectors. Sort the S' eigenvalues ​​in ascending order and extract the c eigenvectors corresponding to the c smallest eigenvalues. The length of each eigenvector is S'.

[0043] In the specific implementation, c feature vectors are arranged in rows to form a feature vector matrix. K-means clustering is performed on each row of the feature vector matrix to obtain the partitioning results of multiple sub-clusters. c column vectors of length S' are concatenated horizontally to form an S' x c feature vector matrix. The e-th row of the feature vector matrix represents the c-dimensional feature representation of the e-th node. The k-means clustering algorithm is executed on all S' row vectors of the feature vector matrix, and the number of clusters is set to c. The k-means clustering algorithm randomly selects c row vectors from the S' row vectors of the feature vector matrix as initial cluster centers. The Euclidean distance between each row vector and each cluster center is calculated, and each row vector is assigned to the cluster containing the nearest cluster center. The cluster center of each cluster is updated to the arithmetic mean of all row vectors within that cluster. This assignment and update step is repeated until the change in cluster centers is less than a preset convergence threshold or a preset maximum number of iterations is reached. The preset convergence threshold is set to 1 x 10^-4, and the preset maximum number of iterations is set to 100. After clustering is completed, each node is assigned a cluster label, with values ​​ranging from 1 to c. Nodes with the same cluster label belong to the same sub-cluster. The c sub-clusters are the result of dividing multiple sub-clusters with similar dynamic behavior characteristics in the material particle group.

[0044] Example 4: In specific implementation, please refer to Figure 5 The orbital distribution topology is unfolded into a one-dimensional time series, containing topological feature values ​​arranged in the acquisition sequence. The orbital distribution topology is an S-row, k-column matrix, where the s-th row represents the coordinates of the phase space point corresponding to the s-th time frame in the low-dimensional manifold. The value of s ranges from 1 to S, where S is the total number of phase space points. For each phase space point in the orbital distribution topology, a topological feature value is calculated. The topological feature value of the s-th phase space point is taken as the difference between the maximum and minimum values ​​of all k coordinate components in the k-dimensional coordinate vector of the s-th phase space point, i.e., the range of the k coordinate components. When all k coordinate components are equal, the topological feature value of the s-th phase space point is 0. The S topological feature values ​​are arranged in the chronological order of the time frames corresponding to the phase space points, with the time frame order consistent with the acquisition timestamp order of the high-dimensional phase space vector sequence, forming a one-dimensional time series of length S. The s-th element in the one-dimensional time series is the topological feature value of the s-th phase space point.

[0045] In practice, a trajectory matrix is ​​constructed for the one-dimensional time series. The row elements of the trajectory matrix consist of subsequences truncated by a sliding window. A sliding window length parameter U is set, where U is an integer greater than 1 and less than S divided by 2. The sliding window length parameter U is determined based on the expected periodic fluctuation scale during material blending. The sliding window length parameter U is set to the integer value obtained by rounding up twice the delay step number corresponding to the first zero-crossing point of the autocorrelation function in the one-dimensional time series to be analyzed. The sliding window starts from the first element of the one-dimensional time series and slides forward with a step size of 1, truncating a continuous subsequence of length U as a row of the trajectory matrix each time. The sliding window slides a total of S-U+1 times, generating a trajectory matrix with S-U+1 rows and U columns. The element value in the g-th row and h-th column of the trajectory matrix is ​​equal to the topological eigenvalue of the (g+h-1)-th element in the one-dimensional time series, where g ranges from 1 to S-U+1 and h ranges from 1 to U. The trajectory matrix is ​​a Hankel matrix with (S-U+1) rows and U columns, and the elements on the anti-diagonal of the trajectory matrix have the same value.

[0046] In practice, singular value decomposition is performed on the trajectory matrix to generate a sequence of singular values ​​and corresponding left and right singular vectors. The trajectory matrix is ​​denoted as X, with S-U+1 rows and U columns. The product of the transpose of the trajectory matrix X and the trajectory matrix X'' is calculated to obtain a U-row, U-column square matrix. The opposing team Perform eigenvalue decomposition to obtain U eigenvalues ​​and their corresponding U orthogonal eigenvectors, where all U eigenvalues ​​are non-negative real numbers. Sort the U eigenvalues ​​in descending order, and denote the j-th eigenvalue as... The value of j is an integer ranging from 1 to U. Calculate the j-th singular value after sorting. Singularity Equal to eigenvalue The arithmetic square root, that is All U singular values ​​are arranged in order of j from 1 to U to form a singular value sequence. The right singular vector is taken as a square matrix. The eigenvectors, the j-th right singular vector after sorting. It is related to eigenvalues The corresponding U-dimensional column vector. For each singular value... For j values ​​greater than 0, calculate the j-th left singular vector. Left singular vector Multiply the trajectory matrix X by the right singular vector Divide by the singular value To obtain, that is Left singular vector It is a (S-U+1) dimensional column vector. All U singular values ​​and their corresponding U left singular vectors and U right singular vectors constitute the singular value decomposition result.

[0047] In practical implementation, the trend component and detail component are reconstructed based on the left and right singular vectors. The trend component reflects the uniformity baseline, and the detail component reflects local fluctuations. The U singular values ​​in the singular value sequence are grouped, with the first singular value being grouped into the next group. The singular value is divided into trend singular value groups, and the () The singular values ​​from the +1th singular value to the Uth singular value are divided into detail singular value groups. This is the cutoff point for the trend singular value group. The value is the singularity index at which the ratio of the first singular value to the second singular value first exceeds the preset trend cutoff ratio, which is set to 10. The trend component matrix is ​​reconstructed using the singular values ​​corresponding to the trend singularity group and the left and right singular vectors. The trend component matrix is ​​a (S-U+1) row, U column matrix. The element value in the g-th row and h-th column of the trend component matrix is ​​given by the following formula:

[0048]

[0049] in, This represents the value of the element in the g-th row and h-th column of the trend component matrix. Indicates the cutoff point of the trend singular value group. This represents the m-th singular value after sorting, where m ranges from 1 to... integers, Represents the m-th left singular vector The g-th component in the equation, where g takes the value of an integer from 1 to S-U+1. Represents the m-th right singular vector The h-th component in the trend component matrix is ​​defined, where h is an integer ranging from 1 to U. The detail component matrix is ​​reconstructed using the singular values ​​corresponding to the detail singular value group and the left and right singular vectors. The detail component matrix is ​​a (S-U+1) row, U column matrix. The value of the element in the g-th row and h-th column of the detail component matrix is ​​calculated by summing all singular values ​​belonging to the detail singular value group using the same product method as in the trend component matrix reconstruction formula, with the summation range from m to... +1 to U. The reconstructed trend component matrix and detail component matrix are subjected to diagonal averaging. The diagonal averaging process is applied to both the trend component matrix and the detail component matrix. The diagonal averaging process transforms the (S-U+1) row, U column matrix into a one-dimensional sequence of length S. The transformation is achieved by calculating the arithmetic mean of all elements in the matrix that satisfy the condition that the element in the g-th row and h-th column has (g+h-1) equal to the same value t. This arithmetic mean is used as the t-th element value of the transformed one-dimensional sequence, where t is an integer ranging from 1 to S. After diagonal averaging, a one-dimensional sequence of length S for the trend component and a one-dimensional sequence of length S for the detail component are obtained. The t-th element value in the trend component one-dimensional sequence represents the uniformity baseline value at the t-th position in the acquisition time series, and the t-th element value in the detail component one-dimensional sequence represents the local fluctuation value at the t-th position in the acquisition time series.

[0050] In practice, spatiotemporal coordinate points where the amplitude of the detail components exceeds a threshold are marked as regions of local component concentration anomalies. The one-dimensional sequence of the detail components contains S local fluctuation values. The standard deviation of these S local fluctuation values ​​is calculated, and the standard deviation of the detail components is denoted as... , The threshold is equal to the square root of the sum of the squares of the differences between each of the S local fluctuation values ​​and the mean of those S local fluctuation values, divided by S. Set the standard deviation of detail components Multiply by a preset threshold factor ,Right now Preset threshold factor The value is 3, which is the preset threshold multiple factor. The value of 3 is based on the three-standard-deviation principle of the normal distribution. When the local fluctuation value follows a normal distribution, the probability of exceeding three standard deviations is less than 0.3%, which can be considered a significant deviation from the normal fluctuation range. Iterate through each element in the one-dimensional sequence of detail components. If the absolute value of the t-th element in the one-dimensional sequence of detail components is greater than the threshold value... If the time index t is not an anomaly, the spatiotemporal coordinate point corresponding to the t-th element value is marked as an anomaly. The one-dimensional time index t is mapped back to the spatial position during the material blending process. The spatial position is determined by the product of the acquisition timestamp and the velocity of the material to be tested. All continuous spatiotemporal coordinate points marked as anomalies are aggregated into a set of one or more continuous regions. The aggregation method is to group the anomalous spatiotemporal coordinate points with adjacent time indices t and adjacent spatial positions into the same region. Several local component concentration anomaly regions are obtained through regional connectivity analysis. Each local component concentration anomaly region consists of a set of spatiotemporal coordinate points.

[0051] Example 5: In specific implementation, please refer to Figure 6 For each local component concentration anomaly region, binarization is performed to generate a binary anomaly region mask. Based on the spatiotemporal coordinates contained in the local component concentration anomaly region, the spatial span and corresponding temporal span of the local component concentration anomaly region in the material conveying direction are determined. The spatial span is equally spaced into a spatial grid, with the grid resolution set to one-tenth of the diameter of the monochromatic laser beam. The temporal span is divided into a temporal grid according to the acquisition timestamp interval of the multi-angle laser scattering signal set. A two-dimensional matrix is ​​constructed with the number of rows equal to the number of spatial grid divisions and the number of columns equal to the number of temporal grid divisions. The row index of the two-dimensional matrix corresponds to the spatial position, and the column index corresponds to the temporal position. All spatiotemporal coordinates in the local component concentration anomaly region are traversed. For each spatiotemporal coordinate point, the row index of the spatial grid and the column index of the temporal grid in which the spatiotemporal coordinate point falls are determined, and the element values ​​at the corresponding row and column indices in the two-dimensional matrix are set to 1. After traversing all spatiotemporal coordinates, all element values ​​in the two-dimensional matrix that are not set to 1 are set to 0. A two-dimensional matrix whose element values ​​are only 0 or 1 is a binary outlier mask.

[0052] In practice, the binary anomaly region mask is covered with grids of different scales, and the minimum number of grids covering the binary anomaly region at each grid scale is recorded. A set of grid scale values ​​is set. Q represents the total number of grid scales, and Q is an integer greater than or equal to 6. Grid scale value This represents the number of spatial grid cells contained within the side length of the square covering grid, where z is an integer ranging from 1 to Q. The first grid scale... Set to 1. For each integer z from 1 to Q-1, set the grid scale. Set as Multiply by the scaling factor The integer value rounded up, the scaling factor The value is set to 2, ensuring that the grid scale value sequence forms a geometric sequence with a common ratio of 2. The grid scale value sequence continues to be generated up to... The process stops when the distance to the mask in the binary anomaly region exceeds half the maximum side length. The final value of Q is the total number of mesh scales that satisfy the stopping condition. For each mesh scale... The binary outlier mask is divided into non-overlapping regions of size . OK A square grid of columns, when the number of rows or columns of the binary outlier mask is not... When the value is an integer multiple, it is insufficient at the mask boundary of the binary outlier region. row or The remaining portion of the column is treated as a single grid block. Iterate through all grid blocks; if a grid block contains at least one matrix element with a value of 1, then that grid block is counted as covering the binary outlier region. Scale the grid. The total number of grids covering the binary anomaly region obtained from the following statistics is recorded as the minimum number of grids. Perform grid coverage and statistical operations on all Q grid scales to obtain Q sets of data pairs, each set containing one grid scale. and the corresponding minimum number of grid cells .

[0053] In practice, a linear regression is performed with the logarithm of the grid scale as the independent variable and the logarithm of the minimum grid number as the dependent variable. The slope of the fitted line is then calculated as the fractal dimension. The natural logarithm is taken for Q sets of data, and the independent variable values ​​are calculated. and dependent variable values ,in, This represents the logarithm with base e. It is calculated using Q data points. On the scatter plot, a straight line is fitted using the least squares method, and the equation of the line is in the form of: Let D be the slope of the fitted line and C be the intercept of the fitted line. The slope D of the fitted line is calculated using the following formula:

[0054]

[0055] Where D represents the fractal dimension of the local component concentration anomaly region, and Q represents the total number of grid scales. Let represent the natural logarithm of the z-th grid scale. This represents the natural logarithm of the minimum number of grid cells required to cover the binary anomaly region at the z-th grid scale. Represents all Q items The sum of values, Represents all Q items The sum of values, Represents all Q items and The sum of products Represents all Q items The sum of squares. The calculated D value is used as the fractal dimension of the local component concentration anomaly region. If the calculated D value is negative, the fractal dimension is set to 0.

[0056] In practice, the fractal dimension is compared with the fractal dimension of the standard homogeneous material to generate a dimensionless exponent. The fractal dimension of the standard homogeneous material is denoted as... The fractal dimension of a standard homogeneous material The method for obtaining the sample is as follows: Each component raw material is weighed using the same component types and proportions as the material to be tested. Under laboratory conditions, each component raw material is placed in a mixing container and thoroughly mixed at a preset stirring speed and duration to obtain a standard homogeneous material. Using the same monochromatic laser source parameters, ring array laser receiver configuration, acquisition timing, and signal processing methods as the material to be tested, the standard homogeneous material undergoes multi-angle laser scattering signal acquisition, multispectral fusion feature map generation, orbital distribution topology extraction, local component concentration anomaly region location, and fractal dimension calculation. The calculated fractal dimension is used as the fractal dimension of the standard homogeneous material. And the fractal dimension of standard homogeneous materials Stored in the data processing system. For each local component concentration anomaly region of the test material, the fractal dimension D of the local component concentration anomaly region is divided by the fractal dimension of the standard homogeneous material. To obtain the ratio Where U represents the dimensionless index of the overall uniformity of the material blending. The closer the dimensionless index U is to 1, the closer the degree of uniformity of the material being tested is to the standard uniform state.

[0057] The above description is merely a specific embodiment of this application, but the scope of protection of this application is not limited thereto. Any changes or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in this application should be included within the scope of protection of this application.

Claims

1. A method for real-time monitoring of material mixing uniformity based on laser scattering, characterized in that, include: Collect the collection of laser scattering signals from multiple angles after the laser beam penetrates the material under test; The multi-angle laser scattering signal set is preprocessed to generate a multispectral fusion feature map, including: performing time-series segmentation on the laser scattering signal at each angle to obtain scattering intensity subsequences with equal time intervals; extracting statistical features from the scattering intensity subsequences, the statistical features including peak intensity, full width at half maximum (FWHM), and skewness coefficient; encoding the statistical features into a two-dimensional matrix in spatial angular order to obtain an initial feature map; and performing three-channel pseudo-color mapping on the initial feature map, where the first channel corresponds to peak intensity, the second channel corresponds to FWHM, and the third channel corresponds to skewness coefficient, to generate the multispectral fusion feature map. The multispectral fusion feature map is input into a chaotic attractor reconstructor to extract the orbital distribution topology of the particle group. This includes: performing pixel-level normalization on the multispectral fusion feature map to generate a normalized feature map; arranging the normalized feature maps in chronological order to form a high-dimensional phase space vector sequence; performing a delayed embedding transformation on the high-dimensional phase space vector sequence to reconstruct the phase space trajectory of the nonlinear dynamic system. The delayed embedding transformation maps the high-dimensional phase space vector sequence to a phase space trajectory matrix, where each row of the phase space trajectory matrix represents a phase space point, and the total number of phase space points is S; extracting the neighborhood connection matrix of the phase space trajectory, which records the connection relationship between each phase space point and its nearest neighbor; setting a low-dimensional manifold dimension parameter k, where k is a positive integer greater than 1 and less than S, representing the number of spatial dimensions after projection; and for each phase space point a in the neighborhood connection matrix, finding the nearest neighbor of phase space point a based on the connection relationship in the a-th row of the neighborhood connection matrix. Then find the nearest neighbor point in the phase space trajectory matrix. The k-1 phase space points with the smallest Euclidean distance, and the phase space point a and its nearest neighbor. Together, they form a set of k neighboring points of phase space point a; using the k neighboring points of phase space point a, a linear reconstruction of phase space point a is performed, and the linear reconstruction weight coefficients are solved; the linear reconstruction weight coefficients of the S phase space points are combined into an S-row S-column weight matrix, where the elements corresponding to non-neighboring points in the weight matrix have a value of 0; based on the weight matrix, a sparse symmetric matrix M is constructed, which is an S-row S-column square matrix, calculated as M = (identity matrix - weight matrix). T (Identity matrix - weight matrix); Perform eigenvalue decomposition on the sparse symmetric matrix M, and extract k eigenvectors corresponding to the 2nd to k+1th eigenvalues ​​of the sparse symmetric matrix from smallest to largest. The length of each eigenvector is S; Arrange the k eigenvectors into an S-row k-column matrix, which is the low-dimensional manifold projection result. Each row of the low-dimensional manifold projection result is a k-dimensional vector, representing the coordinates of a phase space point in the low-dimensional manifold. The set of all S k-dimensional vectors is the orbital distribution topology. Singular spectrum analysis was performed on the orbital distribution topology to locate local component concentration anomalies during the material blending process; Calculating the statistical self-similarity parameters of the local component concentration anomaly regions and constructing a dimensionless index representing the overall blending uniformity of the material includes: binarizing each local component concentration anomaly region to generate a binary anomaly region mask; covering the binary anomaly region mask with grids of different scales and recording the minimum number of grids covering the binary anomaly region at each grid scale; performing linear regression with the logarithm of the grid scale as the independent variable and the logarithm of the minimum number of grids as the dependent variable, and calculating the slope of the fitted line as the fractal dimension; and calculating the ratio of the fractal dimension to the fractal dimension of the standard homogeneous material to generate the dimensionless index.

2. The method for real-time monitoring of material blending uniformity based on laser scattering according to claim 1, characterized in that, The singular spectrum analysis of the orbital distribution topology to locate local component concentration anomalies during material blending includes: The orbital distribution topology is unfolded into a one-dimensional time series, which contains topological feature values ​​arranged in the order of acquisition time. A trajectory matrix is ​​constructed for the one-dimensional time series, and the row elements of the trajectory matrix are composed of subsequences truncated by a sliding window; Perform singular value decomposition on the trajectory matrix to generate a sequence of singular values ​​and corresponding left and right singular vectors; The trend component and detail component are reconstructed based on the left and right singular vectors, wherein the trend component reflects the uniformity baseline and the detail component reflects local fluctuations. The spatiotemporal coordinate points where the amplitude of the detailed component exceeds the threshold are marked as the local component concentration anomaly region.

3. The method for real-time monitoring of material blending uniformity based on laser scattering according to claim 1, characterized in that, The collection of multi-angle laser scattering signals after the laser beam penetrates the material under test includes: A ring-shaped laser receiver is provided, which includes multiple photodetectors distributed at equal angular intervals. A monochromatic laser source is excited and incident on the test material in a vertical direction, causing the test material to generate scattered light; The ring array laser receiver is synchronously triggered, so that each photodetector simultaneously collects the scattering intensity data at the corresponding angle; The scattering intensity data is encapsulated according to angle index and acquisition timestamp to generate the multi-angle laser scattering signal set.

4. The method for real-time monitoring of material blending uniformity based on laser scattering according to claim 1, characterized in that, Before inputting the multispectral fused feature map into the chaotic attractor reconstructor, the method further includes: The multispectral fusion feature map is subjected to blind source separation processing to remove the common-mode noise component introduced by ambient stray light; Median filtering is performed on abnormal pixels in the multispectral fusion feature map to repair signal loss areas caused by dust or bubbles. The multispectral fusion feature map is spatially downsampled to obtain a uniform-sized input feature map.

5. The method for real-time monitoring of material blending uniformity based on laser scattering according to claim 1, characterized in that, After inputting the multispectral fused feature map into the chaotic attractor reconstructor to extract the orbital distribution topology of the material particle group, the process further includes: Calculate the global Lyapunov exponent of the orbital distribution topology, whereby the global Lyapunov exponent characterizes the degree of chaos in the material mixing system; By deleting isolated nodes from the neighborhood connection matrix of the track distribution topology, a simplified topology is obtained. Spectral clustering analysis was performed on the simplified topology to divide the material particle group into multiple sub-clusters with similar dynamic behavior characteristics.

6. The method for real-time monitoring of material blending uniformity based on laser scattering according to claim 4, characterized in that, The blind source separation process uses an independent component analysis algorithm, with the number of independent components set to 3.

Images