Automatic injection molding product detection method and system based on machine vision

Abstract

The invention relates to the technical field of injection molding product quality automatic detection, and discloses an injection molding product automatic detection method and system based on machine vision, and the method comprises the steps: obtaining a polarization camera image sequence and an event camera data flow which are synchronously triggered, aligning the two kinds of data according to timestamps to generate a multi-modal original data set; fusing the standardized phase change rate and event density based on a Kalman filter, and generating a stress-vibration coupling feature vector containing a stress state and a vibration state through state estimation; analyzing a spatio-temporal evolution mode of the stress-vibration coupling feature vector, identifying the stress-vibration coupling feature vector as a crack propagation precursor signal when the spatial gradient and the time derivative of the feature vector simultaneously exceed preset thresholds, and outputting failure position prediction and residual life evaluation results according to the spatial position and the feature evolution rate of the precursor signal. According to the invention, accurate monitoring and failure prediction of the whole crack propagation process are realized.

Description

Technical Field

[0001] This invention relates to the field of automatic quality inspection technology for injection molded products, and more specifically, to an automatic inspection method and system for injection molded products based on machine vision. Background Technology

[0002] In the field of quality inspection of medical injection molded devices, especially in the cyclic pressure resistance testing of transparent devices such as syringe barrels and infusion set tee connectors, it is necessary to monitor the structural integrity of the devices in real time during dynamic pressurization. However, in the cyclic pressure testing of medical devices, the initiation and propagation of microcracks is a complex dynamic process involving slow stress accumulation and instantaneous stress release. Although traditional polarization imaging technology can detect stress distribution, its frame rate is usually below 100fps, which cannot capture the microsecond-level stress release vibrations generated during the instantaneous propagation of cracks; while using high-speed cameras or vibration sensors alone cannot obtain spatial distribution information of the stress field, resulting in the inability to accurately predict the location and timing of failure. Summary of the Invention

[0003] This invention provides an automatic inspection method and system for injection-molded products based on machine vision, overcoming the limitations of single imaging modal methods in terms of temporal resolution or spatial information completeness in related technologies. Traditional methods cannot simultaneously capture slow stress changes and transient vibration characteristics. This addresses the technical problem of not being able to simultaneously obtain high temporal resolution and complete spatial information in dynamic stress monitoring of medical devices.

[0004] This invention provides an automatic inspection method for injection-molded products based on machine vision, comprising:

[0005] Acquire synchronously triggered polarization camera image sequences and event camera data streams, extract polarization phase information from the polarization camera image sequences, obtain pixel-level brightness change events from the event camera data streams, and align the two types of data by timestamps to generate a multimodal raw dataset;

[0006] The mutual information maximization algorithm is used to spatially register the polarization image coordinate system and the event camera coordinate system. The coordinate transformation matrix is ​​solved by maximizing the mutual information between the edge features of the polarization image and the cumulative event map.

[0007] Calculate the phase difference between adjacent moments of the registered polarization phase data. When the phase difference exceeds the phase change threshold, it is marked as a phase change region. The event density of the corresponding spatiotemporal location of the phase change region is then calculated.

[0008] Based on the normalized phase change rate and event density fused by the Kalman filter, a stress-vibration coupling feature vector containing stress state and vibration state is generated through state estimation;

[0009] The spatiotemporal evolution pattern of stress-vibration coupled eigenvectors is analyzed. When the spatial gradient and time derivative of the eigenvectors both exceed a preset threshold, they are identified as precursor signals for crack propagation. Based on the spatial location and eigenvalue evolution rate of the precursor signals, the failure location prediction and remaining life assessment results are output.

[0010] Furthermore, the extraction of the polarization phase information includes: calculating the Stokes parameters from images at at least four different polarization angles, and calculating the phase value using the four-quadrant arctangent function to ensure the accuracy and continuity of the phase information across the entire quadrant.

[0011] Furthermore, the mutual information maximization algorithm includes: extracting edge features of the polarization image, accumulating events within the time window into a two-dimensional image, performing spatial transformation on the accumulated event map through a parameterized homography transformation model, and optimizing the mutual information objective function using the gradient ascent method;

[0012] The homography matrix satisfies the following constraints: the determinant of the matrix is ​​greater than zero, the scale constraint is between 0.5 and 2.0, the translation constraint is no more than one-quarter of the image size, and the rotation constraint is limited by singular value decomposition to ensure that the ratio of the maximum to the minimum singular value does not exceed 5.

[0013] Furthermore, the calculation of the event density includes: counting the number of events within a set time window within a spatial neighborhood centered on the target location and with a set radius;

[0014] The phase change threshold is determined based on statistical analysis of phase changes under normal pressure conditions and is set as a multiple of the standard deviation of the phase difference. The multiple factor is selected between 2 and 4 according to the instrument type and detection accuracy requirements.

[0015] Furthermore, the state vector of the Kalman filter includes two components: stress state and vibration state, and the observation vector consists of the standardized phase change rate and event density.

[0016] Among them, the standardization preprocessing adopts the min-max standardization method to normalize the phase change rate and event density to the 0 to 1 range, respectively, to eliminate the influence of different physical dimensions on the fusion process.

[0017] Furthermore, before calculating the phase difference between adjacent time points, the algorithm also includes: processing the polarization phase map under each pressure state based on the least squares phase unwrapping algorithm, restoring the continuous phase distribution by minimizing the inconsistency of the phase gradient, and ensuring the spatiotemporal continuity of the phase information.

[0018] Furthermore, after fusion based on Kalman filters, the method also includes: inputting phase change region data and event cluster features into a dual-flow graph neural network; the stress topology flow constructs the phase change region into a graph structure to handle stress transmission relationships; the vibration propagation flow constructs the event cluster into a dynamic graph to handle vibration propagation paths; and the two flows interact with each other through a cross-flow attention mechanism.

[0019] Furthermore, the remaining life assessment is based on the Weibull distribution model. By collecting historical failure data of similar devices, the maximum likelihood estimation method is used to fit the Weibull distribution parameters, and the remaining life is calculated based on the normalized spread rate and the Weibull distribution parameters.

[0020] Furthermore, before analyzing the spatiotemporal evolution pattern, the method also includes: inputting the coupled feature vector sequence into a long short-term memory network, learning the temporal evolution law of stress accumulation and release through the gating unit update mechanism, outputting the feature prediction value of future time, and improving the temporal accuracy of failure prediction.

[0021] This invention provides an automatic inspection system for injection molded products based on machine vision, comprising:

[0022] The multimodal image acquisition module, including a polarization camera and an event camera, is used to simultaneously acquire polarization image sequences and event data streams of injection-molded products during cyclic pressure testing.

[0023] The spatial registration module is used to achieve accurate registration of the coordinate systems of the polarization camera and the event camera through the mutual information maximization algorithm;

[0024] The feature extraction module is used to extract phase abrupt change regions and calculate event density distribution;

[0025] The feature fusion module includes a Kalman filter for fusing multimodal features to generate stress-vibration coupled feature vectors;

[0026] The failure prediction module is used to analyze the spatiotemporal evolution pattern of feature vectors and output the failure location and remaining life prediction results.

[0027] The beneficial effects of this invention are as follows:

[0028] This invention overcomes the limitations of single imaging modalities in temporal resolution or spatial information integrity by simultaneously acquiring image sequences from a polarization camera and data streams from an event camera, and using a mutual information maximization algorithm for spatial registration. By calculating the phase difference between adjacent time points to extract phase abrupt change regions and analyzing the event density distribution at corresponding locations, precise localization of abrupt events during stress accumulation is achieved. A stress-vibration coupling feature vector is generated by fusing phase change rate and vibration energy distribution using a Kalman filter, organically combining the static stress accumulation process with dynamic stress release events, overcoming the problem that traditional methods cannot simultaneously capture slow stress changes and transient vibration characteristics. Therefore, this invention solves the technical problem of not being able to simultaneously obtain high temporal resolution and complete spatial information in dynamic stress monitoring of medical devices, achieving precise monitoring and failure prediction of the entire crack propagation process. Attached Figure Description

[0029] Figure 1 This is a flowchart of an automatic inspection method for injection molded products based on machine vision according to the present invention;

[0030] Figure 2 This is a bi-line graph of the phase change rate and event density time evolution sequence of the present invention;

[0031] Figure 3 This is a scatter plot comparing the coordinates before and after spatial registration according to the present invention;

[0032] Figure 4 This is a biaxial evolution and bi-Y-axis hybrid diagram of the stress-vibration coupling state of the present invention;

[0033] Figure 5 This is the Sankey diagram of the data flow of the detection system of the present invention. Detailed Implementation

[0034] The subject matter described herein will now be discussed with reference to exemplary embodiments. It should be understood that these embodiments are discussed only to enable those skilled in the art to better understand and implement the subject matter described herein, and changes may be made to the function and arrangement of the elements discussed without departing from the scope of this specification. Various processes or components may be omitted, substituted, or added as needed in the examples. Furthermore, some features described in the examples may be combined in other examples.

[0035] At least one embodiment of the present invention discloses an automatic inspection method for injection-molded products based on machine vision, such as... Figure 1 As shown, it includes the following steps:

[0036] Step 1: Obtain the synchronously triggered polarization camera image sequence and event camera data stream to generate a multimodal raw dataset;

[0037] Image sequences with different polarization angles are acquired from a polarization camera, and polarization phase information is extracted. ,in For image coordinates, For time;

[0038] Simultaneously acquire event data streams from the event camera. ,in For event space coordinates, For event timestamps, Polarity is determined. The two datasets are aligned by timestamps to generate the original multimodal dataset.

[0039] It should be noted that polarization phase information The extraction was achieved by calculating the Stokes parameters from images at at least four different polarization angles (0°, 45°, 90°, 135°).

[0040]

[0041]

[0042] in These are the light intensity values ​​corresponding to the polarization angles.

[0043] Furthermore, to ensure the accuracy and continuity of phase information, the four-quadrant arctangent function atan2 is used for phase calculation:

[0044]

[0045] The output range of the four-quadrant arctangent function atan2 is: It can correctly handle phase angles in all quadrants, avoiding the limitation of the standard arctan function which only covers phase angles in all quadrants. Limitations of scope.

[0046] Step 2: Use the mutual information maximization algorithm to spatially register the polarization image coordinate system and the event camera coordinate system, and calculate the coordinate transformation matrix;

[0047] The edge feature map of the polarization image is registered with the cumulative event map of the event camera within the same time window, and the mutual information between the two is maximized. To solve the coordinate transformation matrix :

[0048]

[0049] in For polarization image edge features, For cumulative event graphs, This is the transformation function.

[0050] Furthermore, the complete constraint condition for the mutual information maximization optimization problem is: the optimization variable is the homography matrix. The constraints include: matrix determinant constraints. Ensure the directionality of the transformation; scale constraints To prevent overscaling, Frobenius norm; translation constraint and The translation amount is limited to no more than 1 / 4 of the image size, where and These represent the image width and height, respectively; rotation constraints are decomposed using SVD. ,Require ,in and These are the maximum and minimum singular values, respectively, to prevent excessive stretching and deformation; the time consistency constraint requires the transformation matrix to be within the time window. The internal stability is maintained, that is ,in Tolerance for time variations.

[0051] The aforementioned mutual information maximization algorithm includes the following processing steps:

[0052] First, extract the edge features of the polarization image. Image gradients are calculated using the Sobel operator. Take the gradient magnitude As edge strength;

[0053] Then the time window The events within are accumulated into a two-dimensional image:

[0054]

[0055] Next, the model is transformed using parameterization. Perform a spatial transformation on the cumulative event graph, where for The homography matrix, the transformation function is defined as ,

[0056] That is, for the input image Each pixel coordinate Applying homography transformation The transformed image is obtained through bilinear interpolation.

[0057] Finally, the gradient ascent method is used to optimize the mutual information objective function.

[0058] Mutual information is defined as ,in For joint probability distribution, and This represents a marginal probability distribution.

[0059] Furthermore, the specific implementation steps of the gradient ascent method for optimizing mutual information are as follows: Initialize the homography matrix. (3×3 identity matrix), set the learning rate and convergence threshold ; in the In the next iteration, the current transformation is calculated. Mutual information value ,

[0060] Then calculate the gradient of mutual information with respect to the transformation parameters. The partial derivatives are calculated using numerical differentiation methods: ,in For minute perturbations, For the first A matrix where each position is 1 and all other positions are 0;

[0061] Then update the transformation matrix. Repeat the iterations until the convergence condition is met. Or reach the maximum number of iterations .

[0062] Step 3: Calculate the phase difference between adjacent time points and extract the phase abrupt change region;

[0063] Calculate the time difference for the registered polarization phase data:

[0064] ,in The sampling interval. When At that time, the position The regions marked as phase transition regions, in which This is the phase abrupt change threshold. Simultaneously, the event density at the spatiotemporal location corresponding to the phase abrupt change region is calculated. .

[0065] Furthermore, the sampling interval The settings should take into account both the frame rate of the polarization camera and the expected rate of stress change to be detected: for cyclic withstand voltage testing, the stress change period is typically... s, to capture sufficient time detail, set ,in For the pressurization cycle, the typical value range is: ms; at the same time, it is necessary to ensure Not less than the minimum frame interval of a polarization camera, i.e. ,in This refers to the frame rate of the polarization camera.

[0066] Furthermore, the phase abrupt change threshold The method for determining it is as follows:

[0067] Based on statistical analysis of phase changes under normal pressure conditions, the standard deviation of the phase difference is calculated. Set the threshold to ,in As a multiple factor, its typical value range is... The specific values ​​are determined based on the type of instrument and the required testing accuracy. For syringe barrel testing, the recommended values ​​are... Recommended for testing infusion set T-connectors .

[0068] It should be noted that event density The calculation is performed by using... Centered on, with radius Within the neighborhood of, the statistical time window The number of events within is obtained as follows:

[0069]

[0070] in This is an indicator function.

[0071] Furthermore, time window parameters and spatial neighborhood radius The setting principle is: time window The frequency should be determined based on the expected vibration frequency to be detected. For stress release events of medical devices, the typical frequency range is [insert range here]. Hz, therefore set ,in The maximum detection frequency typically ranges from [value range missing]. ms; spatial neighborhood radius The value should be set according to the image resolution and the expected defect size, and is usually set to [value missing]. Pixels ensure that enough spatial neighborhood information is captured while avoiding over-smoothing.

[0072] The processing steps of the aforementioned phase change detection and event density coupling algorithm include:

[0073] First, for each pixel position Calculate the phase-time derivative Then, a binary mutation marker is generated based on a threshold. ;

[0074] Next, morphological processing is performed on the marked regions, and neighboring mutation points are merged through connected component analysis to form mutation regions. Finally, the event response intensity of each mutation region is calculated. ,in The area represents the region. The phase abrupt change detection and event density coupled algorithm outputs a set of abrupt change regions. and the corresponding event response intensity .

[0075] Furthermore, the specific implementation methods for morphological processing and connected component analysis are as follows: First, the binary mutation markers are processed... Perform morphological closing operations to fill small gaps, using circular structural elements The closing operation is defined as follows: ,in This represents the expansion operation. This represents the erosion operation;

[0076] Then, 8-neighbor connectivity is used for connected component labeling for each pixel. Check its 8 neighboring pixels Pixels with a value of 1 that are interconnected are grouped into the same connected component; then, connected components with excessively small areas are filtered out, and a minimum area threshold is set. Pixels, the area to be deleted is smaller than Connected components are used to remove noise; finally, for each valid connected component... Calculate its geometric characteristics:

[0077] centroid coordinates ,area Aspect Ratio of Bounding Box .

[0078] Step 4: Calculate the stress-vibration coupling eigenvector by fusing the phase change rate and vibration energy distribution using a Kalman filter;

[0079] Phase change rate and event density As an observation, state estimation is performed using a Kalman filter. Since the phase change rate (unit: radians / time) and event density (unit: number of events / time) have different physical dimensions and numerical ranges, the observation data needs to be standardized preprocessed: the phase change rate and event density are normalized to... The interval eliminates the influence of dimensional differences on the fusion process.

[0080] Furthermore, the specific formula for the standardized preprocessing is as follows: The phase change rate is standardized using a minimum-maximum normalization method. ,

[0081] in and These represent the maximum and minimum values ​​of the phase change rate during the observation period, respectively; the same standardization method is used for event density. ,in and These represent the maximum and minimum event densities within the observation period, respectively. Minimum-maximum normalization ensures that both types of observations are mapped to... Within this interval, it facilitates effective fusion by the Kalman filter. The Kalman filter state estimation formula is:

[0082]

[0083] Where the state vector Includes stress state and vibration state Observation vector For standardized observations,

[0084] Here is the state transition matrix. For the observation matrix, Kalman gain. Output stress-vibration coupling eigenvector. .

[0085] Furthermore, the complete definition and constraints of the Kalman filter parameters are as follows:

[0086] State vector constraints and Ensure that the state estimate is within the standardized range; observation period length constraint. ,in To maximize the number of time steps, ensure that sufficient periodic data is included;

[0087] Time window constraints of standardized parameters ,in and To standardize the start and end times of the computation and ensure the stability of the statistics; the stability constraint of the state transition matrix requires eigenvalues To ensure system stability; numerical stability constraints of Kalman gain. ,in For condition numbers, This is a condition number threshold to prevent numerical instability during matrix inversion.

[0088] Prediction error covariance matrix update Posterior covariance matrix update Initial covariance matrix .

[0089] The input to the aforementioned Kalman filter includes: the sequence of observation vectors. , where each observation vector It includes two components: the standardized phase change rate and the event density;

[0090] System model parameters include the state transition matrix. Observation matrix Process noise covariance matrix and observation noise covariance matrix The output of the Kalman filter is the state estimate at each time step. That is, the optimal estimates of the stress state and vibration state. Kalman gain is achieved through... Calculation, where Let be the prior error covariance matrix.

[0091] Furthermore, the process noise covariance matrix and observation noise covariance matrix The method for determining the process noise covariance matrix is ​​as follows: The uncertainty of the system state is reflected by analyzing the statistical changes in the time series of stress-vibration eigenvectors under normal conditions, specifically calculated as follows: ,in Let Variance be the variance of the stress state variation. The variance of the vibration state change is obtained through statistical analysis of the state differences at consecutive time points in historical data, with a typical range of values. , ; Observation noise covariance matrix The noise level of the measurement system is reflected by repeatedly measuring the phase change rate and event density under the same state, and is calculated as follows: ,in To measure the noise variance of the phase change rate, The noise variance is measured for event density, with a typical range of values. , .

[0092] In this embodiment of the application, in order to improve the continuity of phase information, the following steps are included before step 3:

[0093] Step 3a: Process the polarization phase maps under various pressure states using a phase unwrapping algorithm to generate continuous phase distribution data. The least squares phase unwrapping algorithm is used to recover the continuous phase by minimizing the inconsistency of the phase gradient.

[0094]

[0095] in For phase wrapping function, the phase difference is limited to Within the scope, specifically defined as ,in This indicates rounding to the nearest integer.

[0096] Furthermore, the complete constraints and time dimension of the phase unwrapping optimization problem are as follows:

[0097] The optimization variable is a continuous phase distribution. ,in Image size, The time series length is given; the spatial continuity constraint requires that the phase gradient between adjacent pixels satisfy the following conditions. ,in Spatial gradient threshold;

[0098] The time continuity constraint requires that the phase change between adjacent time points satisfy:

[0099] ,in The threshold is the time-varying threshold; the boundary constraint requires that the unwrapped phase remain consistent with the original wrapped phase at the image boundary, i.e., at the boundary. Above ;

[0100] The time-domain optimization objective function is expanded to:

[0101]

[0102]

[0103] in For a moment The weighting coefficients, This is the time regularization coefficient, with typical values... .

[0104] The input to the aforementioned least-squares phase unwrapping algorithm is the wrapped phase map. The entangled phase diagram is calculated using the atan2 function in step 1. The algorithm outputs the unentangled continuous phase distribution. Its range is no longer limited to Within the range. The least squares phase unwrapping algorithm constructs a discrete Poisson equation. The Fast Fourier Transform is used to solve the problem, ensuring the continuity of the phase gradient.

[0105] Furthermore, to ensure consistency between phase unwrapping and subsequent processing steps, the continuous phase distribution after unwrapping... The continuity constraint must be met: the phase difference between adjacent pixels should be less than 1 / 3. ,Right now ,and This ensures the time difference calculation in subsequent step 3. It can accurately reflect the true phase changes.

[0106] In this embodiment of the application, in order to more accurately analyze the spatial topology and vibration propagation characteristics of the stress field, the following steps are included after step 4:

[0107] Step 5: Input the phase abrupt change region data and event cluster features into the dual-flow graph neural network to process the stress topology and vibration propagation path, respectively. Since the phase abrupt change region features and event cluster features originate from different physical measurement processes and have different dimensions and numerical distributions, the input features need to be standardized preprocessed: all feature vectors are normalized to the same numerical range, and one-hot encoding is performed on features of different classification types. The dual-flow graph neural network includes two branches: stress topology flow and vibration propagation flow. The stress topology flow constructs the phase abrupt change region into a graph structure. , where nodes Indicates a region of mutation, edge This indicates the stress transfer relationship between regions;

[0108] Vibration propagation flow constructs event clusters as dynamic graphs , where nodes Represents the center of the event cluster, edges Indicates the vibration propagation path, The time dimension is considered. The two streams interact with each other through a cross-stream attention mechanism.

[0109] Furthermore, the event cluster feature extraction employs a density clustering method: firstly, the event density distribution at each time step is calculated. Then, the DBSCAN algorithm is used to cluster high-density regions. Specifically, a density threshold is set. ,

[0110] in and These are the mean and standard deviation of the event density, respectively. Pixels with a density exceeding the threshold are selected as the candidate event point set. The candidate point set is clustered using the DBSCAN algorithm, with the neighborhood radius set to 1. (in (Spatial radius for event density calculation), minimum number of points set to The centroid of each cluster serves as the center of the event cluster, and the event cluster feature vector includes: centroid coordinates. Number of events within a cluster Average event intensity and time tags .

[0111] Furthermore, the determination of stress transfer relationships is based on the spatial proximity and phase similarity of the phase abrupt change regions:

[0112] For any two mutation regions and Calculate the Euclidean distance between their centroids:

[0113] ,when At that time, it was assumed that there was a spatial proximity relationship, in which Set it to 1 / 20 of the image diagonal length; simultaneously calculate the average phase difference between the two regions. ,in For the region The average phase value of the inner pixel, when Phase similarity is assumed to exist; only when both spatial proximity and phase similarity conditions are met is a condition set in the adjacency matrix. ,otherwise .

[0114] Furthermore, the construction of vibration propagation paths is based on the spatiotemporal correlation of event clusters: for event clusters at different times... and ,

[0115] When time continuity is satisfied (in Set to the maximum event propagation delay, typically a value of ) and spatial continuity (in The maximum vibration propagation speed is set according to the material properties of the instrument; a typical value is 1 / 10 of the sound velocity in the material.

[0116] Set in the adjacency matrix of the dynamic graph ,in The distance weight parameter is set to 1 / 2 of the average event cluster spacing.

[0117] In the aforementioned dual-flow graph neural network, the input layer of the stress topological flow receives the node feature matrix. and adjacency matrix ,in The stress characteristic dimension; the feature matrix of the input layer receiving nodes of the vibration propagation flow. Adjacency matrix and time encoding matrix ,in For vibration characteristic dimensions, The time-encoded dimension.

[0118] Furthermore, the principle for setting the network layer dimension parameters is: stress feature dimension Determined based on the number of features in the phase abrupt change region, including the region's centroid coordinates. Average phase value Area Shape parameters (Aspect Ratio) and Phase Gradient Magnitude Therefore, set Vibration characteristic dimension Determined based on the number of features of the event cluster, including centroid coordinates. Number of events within a cluster Average event intensity Time tags and cluster shape parameters Therefore, set ;

[0119] Time encoding dimension A positional encoding method is used, which encodes time information using sine and cosine functions, with the encoding dimension set to [value missing]. The specific coding formula is as follows:

[0120] , ,in For time location, For encoding dimension indexes.

[0121] Each stream of the aforementioned dual-stream graph neural network contains multiple graph convolutional layers, the first... The graph convolution operation for a layer is as follows:

[0122]

[0123] in To add self-connected adjacency matrices, This is the corresponding degree matrix. For the first Layer node representation, The weight matrix is ​​a learnable matrix. This is the activation function.

[0124] The aforementioned cross-stream attention mechanism achieves feature interaction by calculating the attention weights between the representations of nodes in the two streams:

[0125]

[0126]

[0127] in and The nodes represent stress flow and vibration flow, respectively.

[0128] Furthermore, the specific implementation steps of the cross-flow attention mechanism are as follows: First, the node representations of stress flow and vibration flow are dimensionally aligned. When the feature dimensions of the two flows are different ( The vibrational flow characteristics are mapped to the stress flow characteristic space through a linear projection layer: ,in Let be the projection weight matrix. The projection bias is used; then, the similarity scores between nodes are calculated, using dot product attention: , of which Scaling is applied to prevent gradient vanishing;

[0129] Next, the similarity scores are softmax normalized to obtain the attention weights:

[0130]

[0131] Finally, the weighted fusion result is calculated. ,

[0132] in To integrate weights and control the degree of influence of cross-stream information.

[0133] The output layer of the aforementioned dual-flow graph neural network contains two branches: stress topological flow outputs the failure probability distribution. Vibration propagation, flow output, and vibration intensity prediction The network output needs to be decoded and transformed to obtain specific executable detection results: the failure probability distribution is obtained through... The operation determines the node index of the highest-risk failure location, and then maps it back to the original image coordinate system to obtain the specific failure location coordinates; vibration intensity prediction is performed using a threshold. Identify areas of abnormal vibration, among which This is the threshold for vibration anomalies.

[0134] Furthermore, vibration anomaly threshold The method for determining it is as follows:

[0135] Based on the vibration intensity distribution under normal conditions in the training data, the value is set to the 99th percentile of the vibration intensity distribution, i.e. This ensures that 99% of normal vibration signals are excluded from the abnormal detection range, with a typical value range of [value missing]. Specific values ​​need to be determined based on historical data statistics for different equipment types. Supervised learning mode is used during training, employing a combined loss function:

[0136]

[0137] in Cross-entropy loss is used for failure location classification. The mean squared error loss is used for vibration intensity regression. For regularization terms, Here are the weight coefficients. The optimization strategy uses the Adam optimizer, with an initial learning rate of 100%. And use an exponential decay strategy.

[0138] Furthermore, the complete time dimension expression of the loss function and the parameter constraints are as follows:

[0139] The total loss function is expanded to time series form:

[0140]

[0141] in The length of the time series. and These are time-varying weighting coefficients used to adapt to the importance of different testing phases; time regularization term:

[0142] Used to ensure the temporal continuity of prediction results;

[0143] The weighting coefficient constraint is and ,in For time regularization weights, typical values ​​are... ;

[0144] The update strategy for time-varying weights is , ,in and The time decay coefficient, and Based on the weights.

[0145] Furthermore, weighting coefficients The determination is made using a grid search method: the weight coefficients are optimized using a grid search on the validation set, and the search range is set to... , , The search steps are 0.1, 0.1, and 0.01 respectively;

[0146] Specifically, the validation set is divided into 5 subsets for cross-validation, and for each weight combination... Calculate the average performance index of 5-fold cross-validation. ,in To balance the importance of classification and regression tasks, the weight combination that maximizes the average performance metric is selected as the final parameter. A typical optimal weight combination is... , , .

[0147] Step 6: Process the event density data using time-frequency analysis algorithms to generate vibration energy spectra and instantaneous frequency distributions. This involves analyzing the event density time series. Perform a short-time Fourier transform:

[0148]

[0149] in For window functions, For frequency.

[0150] Extracting vibrational energy and instantaneous frequency .

[0151] Furthermore, the complete parameter definition and constraints of the short-time Fourier transform are as follows:

[0152] Input signal The domain is ,in Total observation duration; window function Hanning windows Window length Satisfying constraints To ensure sufficient time resolution; frequency range constraints are as follows: ,in The lowest resolvable frequency, For the Nyquist frequency, The sampling interval is; the specific form of the vibration energy integral is as follows: Integral constraints require and The constraints for calculating instantaneous frequency are: ,

[0153] And require ,

[0154] in The energy threshold is used to ensure that the extracted instantaneous frequency has sufficient signal strength.

[0155] The input to the aforementioned short-time Fourier transform is an event density time series. Event density time series represents each spatial location The change in the occurrence rate of events at a location over time.

[0156] Window function Select the Hanning window and set the window length to... The overlap rate is set to 50%. The algorithm output includes: time spectrum. Represents the complex amplitude at each time frequency point; vibrational energy Represents the total energy at each moment; instantaneous frequency. This represents the dominant frequency component at each moment.

[0157] Furthermore, the short-time Fourier transform window length The setting is based on the principle of balancing time and frequency resolution: a window length that is too short will result in insufficient frequency resolution, while a window length that is too long will result in insufficient time resolution. According to Heisenberg's uncertainty principle, the product of time and frequency resolution has a lower bound. Regarding the frequency characteristics of stress release events in medical devices, the main frequency components are concentrated in... In the Hz range, to obtain sufficient frequency resolution, the window length is set. ,in Hz is the lowest detection frequency. This is the resolution factor, with typical values. correspond s;

[0158] At the same time, it is necessary to ensure that the window length includes a sufficient number of sampling points, i.e. ,in The sampling interval is adjusted when the calculated window length corresponds to an insufficient number of sampling points. .

[0159] Step 5: Based on the spatiotemporal distribution of precursor signals and crack propagation rate, output the failure location prediction and remaining life assessment results;

[0160] Analysis of stress-vibration coupling eigenvectors The spatiotemporal evolution pattern, when the spatial gradient of the feature vector is detected. and time derivative Simultaneously, when the threshold is exceeded, it is identified as a precursor signal for crack propagation. The failure location is predicted based on the spatial location of the precursor signal, and the remaining life assessment result is output based on the characteristic evolution rate and the instrument fatigue life database.

[0161] The aforementioned spatiotemporal evolution pattern analysis algorithm includes the following processing steps: First, calculate the spatial gradient field of the eigenvectors. ,

[0162] Discrete differentiation is performed using the Sobel operator; then the gradient magnitude is calculated. and rate of change over time Next, a precursor signal marker is generated by determining the threshold. ,in and Spatial and temporal thresholds were set respectively; then, cluster analysis was performed on the marked regions to aggregate spatially adjacent precursor signal points into potential failure regions. The cluster analysis uses the DBSCAN algorithm, with the input being a set of precursor signal points. The neighborhood radius is set to The minimum number of points is set to The output is a set of cluster labels, with each cluster corresponding to a potential failure region; finally, the expansion rate of each failure region is calculated. Due to the rate of change over time It is calculated from stress and vibration components of different dimensions. The propagation rate needs to be normalized to eliminate the influence of dimensions. Then, based on the normalized propagation rate and the Weibull distribution parameters in the historical failure database, through... Estimate remaining lifespan, of which The normalized expansion rate, For the initial lifespan, and These are the parameters of the Weibull distribution.

[0163] Furthermore, the Weibull distribution parameters and The method to obtain it is as follows:

[0164] Collect historical failure data of similar medical devices to establish a sample dataset of expansion rate and remaining lifespan. ,in For the first The normalized expansion rate of each sample, This corresponds to the actual remaining lifespan;

[0165] The maximum likelihood estimation method is used to fit the parameters of the Weibull distribution, and the likelihood function is:

[0166] By solving and Obtain parameter estimates; for syringe-type instruments, typical parameter values ​​are... , For infusion set T-connectors and similar devices, the typical parameter values ​​are: , Once the parameters are determined, they are updated periodically based on new failure data, with the update cycle set at every 100 new samples or every 6 months.

[0167] Furthermore, the complete constraints and numerical stability requirements for Weibull distribution parameter estimation are as follows: sample data constraints. To ensure the reliability of statistical estimates; extend rate data constraints. ,in Maximum normalized rate; remaining lifetime data constraints ,in Hour, Hours to ensure the physical validity of the data; parameter search space constraints. and ,

[0168] Typical search range , Log-likelihood function Numerical calculations require preventing underflow. Time cutoff is 700; numerical stability constraints for gradient calculation and ;

[0169] Parameter convergence criterion and ,

[0170] in , The convergence threshold; the maximum number of iterations is limited. To prevent infinite loops.

[0171] Furthermore, spatial gradient threshold and the threshold of the rate of change over time The method for determining it is as follows: By analyzing the statistical distribution of the gradient in the eigenvector space under normal conditions, the value is set as the 95th quantile of the gradient magnitude distribution, with a typical range of values. ; By analyzing the statistical distribution of the time derivative of the eigenvector under normal conditions, and setting it as the 90th quantile of the time rate of change distribution, the typical value range is [value range missing]. In the DBSCAN clustering algorithm, the neighborhood radius... Set to 2-3 times the pixel pitch, minimum number of dots The value is set according to the expected size of the failure area, and is usually between 5 and 10.

[0172] In this embodiment, to improve the accuracy of temporal prediction, the following steps are included before step 5: Step 4a: Input the coupled feature vector sequence into the Long Short-Term Memory (LSTM) network to learn the temporal evolution of stress accumulation and release. The hidden state of the LTM network is updated as follows:

[0173]

[0174]

[0175] in These are the forget gate, input gate, and output gate, respectively. In cellular state, This represents element-wise multiplication. The network outputs the feature predictions for future timeframes. .

[0176] The input layer of the aforementioned Long Short-Term Memory network receives a sequence of stress-vibration coupled feature vectors:

[0177] ,in For each feature vector, the number of historical time steps is [number]. It includes two components: stress state and vibration state. The input sequence is mapped to a high-dimensional feature space through an embedding layer. ,in For the embedded dimension.

[0178] Furthermore, the method for setting the dimensional parameters of the Long Short-Term Memory network is as follows: historical time steps. Based on the time scale of stress accumulation, for the slow stress accumulation process in cyclic withstand voltage testing, the following settings are defined: ,in For the pressurization cycle, a typical value is... ;

[0179] Embedded Dimension Using the expansion factor of the input dimension,

[0180] Set as To ensure sufficient representational power while avoiding overfitting; hidden state dimension Determined based on sequence complexity and available computational resources, using empirical formulas. For typical parameters and Calculations yielded .

[0181] The formula for calculating the gating unit of the aforementioned Long Short-Term Memory (LSTM) network is as follows:

[0182]

[0183]

[0184]

[0185]

[0186] in This is the weight matrix. For bias vectors, It is the sigmoid activation function. For a moment Input.

[0187] The output layer of the aforementioned Long Short-Term Memory network contains a linear mapping layer that stores the hidden state at the last time step. Mapping to the prediction space:

[0188]

[0189] in To output the weight matrix, For output bias, To hide the state dimension, the predicted feature vector output by the Long Short-Term Memory network needs to be decoded into specific stress and vibration state assessments: the first component of the predicted vector is decoded as the stress accumulation level, which is compared with historical stress thresholds to determine whether it is close to the failure critical point; the second component is decoded as vibration activity, which is used to predict the probability of vibration anomalies in the future.

[0190] The aforementioned Long Short-Term Memory (LSTM) network was trained using a sequence-to-sequence supervised learning model, with the loss function being:

[0191]

[0192] in The first term represents the mean squared loss of the prediction error, and the second term represents the weight regularization term. Here is the regularization coefficient. The optimization strategy uses the RMSprop optimizer, with a learning rate set to... Gradient clipping is used to prevent gradient explosion, with the clipping threshold set to 5.0.

[0193] Furthermore, the complete parameter definition and time constraints of the loss function for Long Short-Term Memory (LSTM) networks are as follows:

[0194] Training sample set ,

[0195] in Given the input sequence, For the target predicted value, Historical time steps; prediction time span constraint ,

[0196] Ensure that the predictions are meaningful and do not exceed the model's capabilities; weight regularization term The summation iterates through all Long Short-Term Memory (LSTM) network layers. Frobenius norm; regularization coefficient constraint Typical values ;

[0197] The time consistency constraint of the loss function requires that the change in loss at adjacent prediction times satisfy... ,in Tolerance for loss variation; the specific implementation of gradient pruning is as follows: ,in For the original gradient, The cropping threshold, Let L2 be the gradient norm.

[0198] This embodiment uses a medical device manufacturer's cyclic pressure resistance test on a 10mL disposable syringe barrel made of polypropylene with a wall thickness of 0.8mm as an example. The testing equipment includes a polarization camera (200fps, 1280×1024 resolution), an event camera (1μs temporal resolution), and a pressure control system. The test pressure range is 0-0.5MPa, and the pressurization cycle is 5 seconds. The goal is to detect the initiation location of microcracks and predict the remaining lifespan before the device experiences macroscopic failure.

[0199] Step 1: Implement an example;

[0200] Table 1. Raw data acquired by the polarization camera:

[0201]

[0202] Table 2: Raw data acquired by the event camera:

[0203]

[0204] Step 2: Implement an example;

[0205] Within a 0.1-second time window, the coordinate transformation matrix is ​​calculated using the mutual information maximization algorithm. for:

[0206]

[0207] Table 3. Verification data of spatial registration results;

[0208]

[0209] Figure 3 The effect of the mutual information maximization algorithm on spatial registration is demonstrated. The accuracy of the spatial registration algorithm is verified by comparing the coordinate distribution before and after registration.

[0210] Step 3: Implement an example;

[0211] At the 15th second of the test, a phase abrupt change was detected in the bottom area of ​​the device.

[0212] Table 4. Phase abrupt change detection results:

[0213]

[0214] Figure 2 This demonstrates the trends in phase change rate and event density over time during cyclic withstand voltage testing.

[0215] Step 4: Implement an example;

[0216] The Kalman filter fusion results show the evolution of the stress-vibration coupling eigenvector.

[0217] Table 5. Calculation results of stress-vibration coupling eigenvectors:

[0218]

[0219] Figure 4 The coupled evolution process of stress state and vibration state is shown, corresponding to the Kalman filter fusion result.

[0220] Step 5: Implement an example;

[0221] Based on Weibull distribution parameters , Failure prediction analysis is performed on the detected precursor signals.

[0222] Table 6. Failure Prediction and Remaining Life Assessment Results:

[0223]

[0224] Figure 5 It showcases the complete technical process from raw data acquisition to final failure prediction, and uses data flow visualization to demonstrate the relationships between each processing step and the data transfer process.

[0225] Actual test results verified that the device failed due to a macroscopic crack at 158.3 seconds into the test, and the predicted remaining lifetime at the initial detection time (15.2 seconds) was 127.4 seconds.

[0226] Prediction error is It meets the accuracy requirements for engineering applications.

[0227] The embodiments of the present invention have been described above. However, the embodiments are not limited to the specific implementation methods described above. The specific implementation methods described above are merely illustrative and not restrictive. Those skilled in the art can make more equivalent embodiments under the guidance of the present embodiments, and all of them are within the protection scope of the present embodiments.

Claims

1. An automatic inspection method for injection-molded products based on machine vision, characterized in that, Includes the following steps: Acquire synchronously triggered polarization camera image sequences and event camera data streams, extract polarization phase information from the polarization camera image sequences, obtain pixel-level brightness change events from the event camera data streams, and align the two types of data by timestamps to generate a multimodal raw dataset; The mutual information maximization algorithm is used to spatially register the polarization image coordinate system and the event camera coordinate system. The coordinate transformation matrix is ​​solved by maximizing the mutual information between the edge features of the polarization image and the cumulative event map. Calculate the phase difference between adjacent moments of the registered polarization phase data. When the phase difference exceeds the phase change threshold, it is marked as a phase change region. The event density of the corresponding spatiotemporal location of the phase change region is then calculated. Based on the normalized phase change rate and event density fused by the Kalman filter, a stress-vibration coupling feature vector containing stress state and vibration state is generated through state estimation; The spatiotemporal evolution pattern of stress-vibration coupled eigenvectors is analyzed. When the spatial gradient and time derivative of the eigenvectors both exceed a preset threshold, they are identified as precursor signals for crack propagation. Based on the spatial location and characteristic evolution rate of the precursor signals, the failure location prediction and remaining life assessment results are output.

2. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, The extraction of polarization phase information includes: calculating the Stokes parameters from images at at least four different polarization angles, and calculating the phase value using the four-quadrant arctangent function to ensure the accuracy and continuity of the phase information across the entire quadrant.

3. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, The mutual information maximization algorithm includes: extracting edge features of the polarization image, accumulating events within the time window into a two-dimensional image, performing spatial transformation on the accumulated event map through a parameterized homography transformation model, and optimizing the mutual information objective function using the gradient ascent method. The homography matrix satisfies the following constraints: the determinant of the matrix is ​​greater than zero, the scale constraint is between 0.5 and 2.0, the translation constraint is no more than one-quarter of the image size, and the rotation constraint is limited by singular value decomposition to ensure that the ratio of the maximum to the minimum singular value does not exceed 5.

4. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, The calculation of the event density includes: counting the number of events within a set time window within a spatial neighborhood centered on the target location and with a set radius; The phase change threshold is determined based on statistical analysis of phase changes under normal pressure conditions and is set as a multiple of the standard deviation of the phase difference. The multiple factor is selected between 2 and 4 according to the instrument type and detection accuracy requirements.

5. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, The state vector of the Kalman filter contains two components: stress state and vibration state, and the observation vector is the standardized phase change rate and event density. Among them, the standardization preprocessing adopts the min-max standardization method to normalize the phase change rate and event density to the 0 to 1 range, respectively, to eliminate the influence of different physical dimensions on the fusion process.

6. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, Before calculating the phase difference between adjacent time points, the algorithm also includes: processing the polarization phase map under each pressure state based on the least squares phase unwrapping algorithm, restoring the continuous phase distribution by minimizing the inconsistency of the phase gradient, and ensuring the spatiotemporal continuity of the phase information.

7. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, After Kalman filter-based fusion, the process also includes: inputting phase change region data and event cluster features into a dual-flow graph neural network; the stress topology flow constructs the phase change region into a graph structure to handle stress transmission relationships; the vibration propagation flow constructs the event cluster into a dynamic graph to handle vibration propagation paths; and the two flows interact with each other through a cross-flow attention mechanism.

8. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, The remaining life assessment is based on the Weibull distribution model. By collecting historical failure data of similar devices, the maximum likelihood estimation method is used to fit the Weibull distribution parameters, and the remaining life is calculated based on the normalized spread rate and the Weibull distribution parameters.

9. The automatic inspection method for injection-molded products based on machine vision according to claim 1, characterized in that, Before analyzing the spatiotemporal evolution pattern, the process also includes: inputting the coupled feature vector sequence into a long short-term memory network, learning the temporal evolution law of stress accumulation and release through the gating unit update mechanism, outputting the feature prediction value of future time, and improving the temporal accuracy of failure prediction.

10. A machine vision-based automatic inspection system for injection-molded products, used to execute the machine vision-based automatic inspection method for injection-molded products according to any one of claims 1-9, characterized in that, include: The multimodal image acquisition module, including a polarization camera and an event camera, is used to simultaneously acquire polarization image sequences and event data streams of injection-molded products during cyclic pressure testing. The spatial registration module is used to achieve accurate registration of the coordinate systems of the polarization camera and the event camera through the mutual information maximization algorithm; The feature extraction module is used to extract phase abrupt change regions and calculate event density distribution; The feature fusion module includes a Kalman filter for fusing multimodal features to generate stress-vibration coupled feature vectors; The failure prediction module is used to analyze the spatiotemporal evolution pattern of feature vectors and output the failure location and remaining life prediction results.

Images