Multi-scale retinex and fastica combined infrared thermal image processing method

Abstract

Disclose a multiscale retinex and FASTICA combined infrared thermal image processing method, including the following steps: laser infrared thermal imaging; FASTICA image processing; multiscale Retinex image noise reduction; signal-to-noise ratio calculation.The FASTICA infrared thermal image processing method based on Retinex can enhance the infrared thermal image to display the key information of the image, and the effect of laser infrared nondestructive testing has obvious visual improvement, and other types of image optimization processing can also be carried out by this method.

Description

Technical Field

[0001] This invention relates to infrared thermal imaging processing technology, specifically to a rapid independent component analysis (FASTICA) infrared thermal imaging processing method based on retinex. Background Technology

[0002] Infrared thermography can capture crucial information, such as various material defects, human symptoms, and interference from equipment. However, numerous interfering factors affect the identification and extraction of key information from infrared thermography, resulting in defects such as low resolution, low contrast, and low signal-to-noise ratio. Infrared thermography enhancement algorithms, by correcting infrared images, facilitate further analysis of key information. After extensive optimization research, various algorithms, such as infrared thermography noise reduction and contrast enhancement algorithms, can be combined to comprehensively improve image quality, depending on the specific circumstances. When detecting defects in specific materials such as carbon fiber laminates using infrared thermography, various factors, including the size and depth of defects, the placement of the specimen during the experiment, and interference from the heating laser, can affect the defect detection results. Infrared image enhancement algorithms have become a powerful means to compensate for poor detection results.

[0003] The blurring of defect information in infrared thermal images may be related to their Gaussian nature. To reduce image Gaussianism and improve the independence of key image information, Independent Component Analysis (ICA) is introduced. ICA is a signal processing method for signal separation. It finds the mathematical relationship between the observed and original signals by solving the mixing matrix. The FASTICA algorithm, based on fixed-point iteration, is a fast ICA algorithm that maximizes negative entropy. It uses the maximum entropy principle to approximate negative entropy and optimizes it using a suitable nonlinear function g. Compared to the ordinary ICA algorithm, it converges faster and does not require setting the step size parameter during iteration.

[0004] Defect detection typically requires comparing defective and non-defective regions, and Retinex and its derivative algorithms offer image processing capabilities that preserve color and enhance detail. Retinex is a nonlinear algorithm whose basic idea is that illumination intensity determines the dynamic range of all pixels in the original image, while the inherent properties of the original image are determined by the object's reflectance coefficient; that is, it is assumed that the original image is obtained by multiplying the reflected image and the illuminated image. When processing infrared images, Retinex can reduce image noise while preserving detail. With further research, multi-scale Retinex can address the problem that single-scale Retinex cannot simultaneously achieve dynamic range compression and color consistency. Summary of the Invention

[0005] To address the problems existing in the prior art, this invention provides an infrared thermal imaging processing method combining multi-scale Retinex and FASTICA, comprising the following steps:

[0006] Step 1: Laser Infrared Thermal Imaging

[0007] A laser infrared thermal imaging non-destructive testing system was used to perform laser infrared non-destructive testing on carbon fiber laminates, and a sequence of num frames of infrared thermal images was obtained through laser infrared thermal imaging.

[0008] Step 2: FASTICA Image Processing

[0009] (1) PCA processing of image data;

[0010] The infrared thermal image sequence is whitened using the PCA preprocessing method, as shown in Equation (1). Physical dimensionality reduction is performed on the infrared thermal image sequence acquired by the laser infrared thermal imaging non-destructive testing system, converting the three-dimensional temperature data sequence into a two-dimensional matrix composed of multiple one-dimensional data vectors. The three dimensions are the frame number, row number, and column number, respectively. The values ​​of the elements in the three-dimensional matrix, determined by these three dimensions, are the temperature values ​​corresponding to the pixels at that frame number, row number, and column number. The three-dimensional temperature data sequence is converted into a two-dimensional matrix composed of multiple one-dimensional data vectors:

[0011] transform_image(i,h)=image(i,j,v) (1)

[0012] Where transform_image(.) represents a two-dimensional matrix composed of multiple one-dimensional data vectors of infrared thermograms. For two-dimensional data of the same frame, the second column of the original data is appended to the end of the first column, the third column is appended to the end of the second column, and so on. The two-dimensional data of the same frame becomes a one-dimensional data sequence. In equation (1), i represents the number of image frames, h represents the one-dimensional temperature sequence of a single frame after image dimensionality reduction, and image(.) represents a three-dimensional temperature data sequence. Where i, j, and v represent the number of frames, rows, and columns corresponding to a certain element in the three-dimensional temperature matrix, respectively.

[0013] To facilitate data processing and fitting, mean reduction is required, which involves removing the mean from the aforementioned two-dimensional matrix, transforming it into a zero-mean vector. After mean reduction and covariance calculation, the covariance matrix Cov is obtained.

[0014] Cov=(transform_image(.)-ave)(transform_image(.)-ave) T (2)

[0015] Where Cov represents the covariance matrix, ave represents the mean temperature of transform_image(.) over all frames, and T represents the transpose; the eigenvector of the covariance matrix Cov is the principal component with the largest variance in the PCA output. This principal component has only one dimension, which is the temperature value after removing the frame number.

[0016] The Cov eigenvector of the covariance matrix is ​​obtained, denoted as eig(Cov). The eigenvector eig(Cov) representing the principal component is multiplied by the physical dimension-reduced two-dimensional data matrix transform_image(.) to obtain the PCA-transformed image data PCAimage(.), which is the eigenvector of the covariance matrix Cov mentioned above, that is, the principal component with the largest variance output by PCA.

[0017] PCAimage(h)=eig(Cov)transform_image(.) (3)

[0018] (2) Fixed-point iterative ICA and optimal solution;

[0019] For simplicity, let the image data after the PCA conversion be V:

[0020] V = PCAimage(h) (4)

[0021] (3) Update the weight matrix w by maximizing the negative entropy;

[0022] Introducing the concept of negative entropy, the maximum negative entropy indicates the maximum non-Gaussianity, meaning the data achieves the highest independence; calculating W... T The negative differential entropy of V requires knowing W. T The probability density distribution function of V is complex to solve, so the approximate formula of equation (5) is used to solve for W. T The maximum negative differential entropy of V is shown in equation (5);

[0023] max(E{Vg(w T V)}-E{g′(w T V)}w)→w best (5)

[0024] In the formula, max(.) represents the maximum value of the polynomial; w best represents the weight matrix corresponding to the maximum negative entropy, indicating that the separation of independent components has been completed; g(·) is the nonlinear function for approximation; T represents matrix transpose, and E represents the expected value;

[0025] Based on the optimal solution theory in nonlinear programming, the equality constraint optimization problem is solved through the Kuhn-Tucker condition, and the objective function for optimization is given by equation (5); under the constraint E{(w T V) 2}=||w|| 2 =1, according to the constraint and objective function (5), after passing the Kuhn Tucker condition, we get equation (6), where β is the non-negative coefficient of Kuhn Tucker;

[0026] E{Vg(w T V)}+βw=0 (6)

[0027] The equation to be solved is equation (6), which is then used to obtain equation (7) via Newton's iteration method. Given the objective function (5) and the constraints of the Kuhn-Tucker condition, the weight matrix w is calculated using Newton's iteration method. k The optimization is performed using Newton's iteration method, as shown in equations (7) and (8), and the iteration stops when the ε threshold is met. k Set to w best ;

[0028]

[0029] Among them, w k+1 This represents the iteration weight vector at step k+1;

[0030] (4) When the newly generated w k+1 Satisfy |w k+1 -w k |<ε, where ε is the threshold, then let w best =w k+1 ; Regarding the weight matrix w best Standardization processing, w best ←w best / ||w best ||;||.|| denotes the 2-norm;

[0031] (5) Take the w obtained in step (4) best Multiply by V to obtain the independent image vector after whitening and ICA; restore the dimension of the independent image vector, converting it from a one-dimensional vector to a two-dimensional vector, where the two-dimensional vector represents the number of rows and columns, and the specific two-dimensional vector value is temperature;

[0032] FASTICA uses the Newton iteration method to obtain the weight matrix within the specified error, as shown in Equation (8). It uses the opposite program logic to dimensionality reduction, that is, the method of dimensionality increase, to process images. It extracts the key information of defects by improving data independence. In data analysis, FASTICA can only consider improving data independence. Different nonlinear functions g(·) are required to process different images. The output image data of the FASTICA processing part is set as I(x,y), where I(x,y) is the original image value.

[0033] I(x, y) = w best V (9)

[0034] Step 3: Multi-scale Retinex Image Denoising

[0035] In step 2, after FASTICA from equation (5) to equation (8), the maximum negative entropy weight w is obtained. best After multiplying with the image matrix V (i.e., Equation (9)), the image is further subjected to multi-scale Retinex image denoising processing to obtain an infrared image combining multi-scale Retinex and FASTICA; the detailed texture information of the image is extracted through Retinex theory. The workflow is as follows: an image I(x,y) is decomposed into illumination component L(x,y) and image reflection component r(x,y), where (x,y) is the image pixel position. The illumination component L(x,y) is decomposed into the temporal convolution of the original image value I(x,y) and the low-pass filter function, and the image reflection component r(x,y) can be expanded into the following equation (10);

[0036] r(x,y)=lnI(x,y)-ln[I(x,y)*LPF(x,y)] (10)

[0037] In the formula, I(x,y) is the original image value, and LPF(x,y) represents the low-pass filter function;

[0038] Low-pass filtering of images is usually implemented using a Gaussian kernel, as shown in equation (11);

[0039]

[0040] In equation (11), Gaussian(x,y) is the Gaussian function, λ is the normalization constant, and δ is the wrapping parameter of the Gaussian kernel. The larger the wrapping parameter δ, the larger the low-pass filtering range. The wrapping parameter δ is the scale of Retinex. Single scale means that the wrapping parameter δ has a unique value, and multi-scale means that the wrapping parameter δ has a non-unique value. By weighting r(x,y), based on the strong effectiveness of the Gaussian kernel's influence within 2 or 3 times the wrapping parameter δ, multiple wrapping parameters are set for multi-scale Retinex processing. The low-pass filtering of the image is performed using multiple wrapping parameters in equation (12). For r with different δ values... M The average value of (x, y) is used to improve the low-pass filtering effect, i.e., r is used. M (x, y) replaces r(x, y) in equation (10), i.e., multi-scale, to improve the output illumination components; the multi-scale retinex expansion (12) is obtained by referring to equation (11).

[0041]

[0042] Step 4: Signal-to-noise ratio calculation

[0043] A comparative signal-to-noise ratio (SNR) calculation method is adopted. By calculating the difference between image vectors, the mean square error of all points in the image is calculated to obtain the SNR* value. When the SNR* value is greater than 1, the quality is improved relative to the original image. Let the two image vectors being compared be Y and X, where the function num(Y) is the number of elements in Y.

[0044] A = YX (13)

[0045] MSE=sum(A(:).*A(:)) / num(Y) (14)

[0046] SNR * =10*log 10 (sum(X(:).*X(:)) / MSE / num(Y)) (15)

[0047] Where A represents the multi-channel pixel difference map of the two images, sum(.) represents the summation of all points within the parentheses, A(:) represents all pixels in the pixel difference map, X(:) represents all pixel values ​​in the original image vector, .* represents matrix dot product, and SNR is calculated. * This indicates the signal-to-noise ratio.

[0048] In one embodiment of the present invention, g is y 3 ,tanh(a1y),yexp(-y 2 One of the types of / 2).

[0049] In one specific embodiment of the present invention, ε is set to 0.001.

[0050] In another specific embodiment of the present invention, in step 3, three wrapping parameters are set for the multi-scale Retinex processing.

[0051] In another specific embodiment of the present invention, in step 3, the first surround parameter δ1 is set to 128, the second surround parameter δ2 is set to 256, and the third surround parameter δ3 is set to 512.

[0052] The Retinex and FASTICA algorithms are well-suited for defect detection at the image processing level. Combining these two algorithms and applying them to the field of material defect detection can address the shortcomings of existing detection methods and improve the defect recognition accuracy of carbon fiber composite infrared thermography. This invention's Retinex-based FASTICA infrared thermography processing method can enhance infrared thermography to display key image information, significantly improving the visual performance of laser infrared non-destructive testing. This method can also be used to optimize images from other industries.

[0053] To verify the superiority of the retinex-based FASTICA infrared thermal imaging processing method, a laser infrared thermal imaging non-destructive testing system was constructed using a high-power laser generator, a beam expander, and an infrared thermal imager. The thermal images of the defect distribution in the carbon fiber laminate were processed and verified. Attached Figure Description

[0054] Figure 1 A schematic diagram of a laser infrared thermal imaging non-destructive testing system is shown.

[0055] Figure 2 The diagram shows the dimensions and defect distribution of carbon fiber laminate No. 1;

[0056] Figure 3 This image shows the effect of FASTICA treatment on carbon fiber laminate No. 1.

[0057] Figure 4 The image shows the effect of processing carbon fiber laminate No. 1 using the multi-scale Retinex method;

[0058] Figure 5 The diagram shows the dimensions and defect distribution of carbon fiber laminate No. 2;

[0059] Figure 6 This image shows the effect of FASTICA treatment on carbon fiber laminate No. 2;

[0060] Figure 7 This image shows the effect of combining carbon fiber laminate No. 2 with multi-scale Retinex processing.

[0061] Figure 8 The image shown is frame 300 of the infrared thermogram sequence for carbon fiber laminate No. 1.

[0062] Figure 9 The image shown is frame 300 of the infrared thermogram sequence for carbon fiber laminate No. 2. Detailed Implementation

[0063] This invention performs laser infrared nondestructive testing on carbon fiber laminates to obtain an infrared thermogram sequence of the composite material. Based on the independent component analysis (FASTICA) with the highest negative entropy, the infrared thermograms are subjected to a series of processing steps, including normalization and whitening, to reduce the Gaussianity of the image signal and obtain independent components with significant nondestructive testing effects. The independent component maps are then enhanced using multi-scale Retinex image enhancement. Finally, the detection effect of this combined algorithm is verified through comparison and a special signal-to-noise ratio method. Specific implementation methods are as follows:

[0064] Step 1: Laser Infrared Thermal Imaging

[0065] A laser infrared thermal imaging non-destructive testing system is composed of a high-power laser generator, a beam expander, and an infrared thermal imager (this device is a commonly used infrared imaging device in this invention, used to capture thermal images of carbon fiber laminates). Figure 1 As shown, laser infrared non-destructive testing was performed on the carbon fiber laminate, and a sequence of 500 infrared thermal images was obtained through laser infrared thermal imaging. Figure 2 , Figure 3 These are two screenshots from an infrared thermal image sequence.

[0066] Step 2: FASTICA Image Processing

[0067] Fast Independent Component Analysis (FASTICA) involves three steps: first, data preprocessing, including mean reduction and whitening; second, selecting a non-Gaussianity decision criterion; and finally, a fixed-point iterative algorithm. In conjunction with infrared thermal image processing, this invention employs Principal Component Analysis (PCA) for infrared thermal image preprocessing. The input to the FASTICA algorithm is the infrared thermal image sequence requiring dimensionality reduction. This infrared thermal image sequence is a three-dimensional image temperature data sequence, obtained by converting the aforementioned infrared thermal image sequence into a new data format. Specifically:

[0068] (1) PCA processing of image data;

[0069] The purpose of whitening is to remove redundant information from the input data. The PCA preprocessing method, commonly used in data processing, is used to whiten the infrared thermal image sequence. The specific PCA process is shown in equation (1). Physical dimensionality reduction is performed on the infrared thermal image sequence acquired by the laser infrared thermal imaging experimental equipment. The three-dimensional temperature data sequence (where the three dimensions are the frame number, row number, and column number, respectively, and the values ​​of the elements in the three-dimensional matrix determined by these three dimensions are the temperature values ​​corresponding to the pixel points at those frame, row, and column numbers; therefore, the obtained infrared thermal image sequence is a three-dimensional matrix, and the value of each element in the matrix is ​​the temperature value at that point defined by the frame, row, and column number of that element) is converted into a two-dimensional matrix composed of multiple one-dimensional data vectors. In this invention, the specific method for converting the three-dimensional temperature data sequence into a two-dimensional matrix composed of multiple one-dimensional data vectors is as follows:

[0070] transform_image(i,h)=image(i,j,v) (1)

[0071] Where transform_image(.) represents a two-dimensional matrix composed of multiple one-dimensional data vectors of a specific infrared thermogram. For two-dimensional data of the same frame, the second column of the original data is appended to the end of the first column, the third column is appended to the end of the second column, and so on. The two-dimensional data of the same frame becomes a one-dimensional data sequence. In equation (1), i represents the number of image frames, h represents the one-dimensional temperature sequence of a single frame after image dimensionality reduction, and image(.) represents the three-dimensional temperature data sequence, where i, j, and v represent the number of frames, rows, and columns corresponding to a certain element in the three-dimensional temperature matrix, respectively. transform_image(i,h)=image(i,j,v) is implemented by setting up three nested for loops in the MATLAB program.

[0072] To facilitate data processing and fitting, mean reduction is required, which involves removing the mean from the aforementioned two-dimensional matrix to obtain a zero-mean vector. The methods and specific implementations of mean reduction are well-known to those skilled in the art and will not be elaborated upon further. After mean reduction and covariance calculation, the covariance matrix Cov is obtained.

[0073] Cov=(transform_image(.)-ave)(transform_image(.)-ave) T (2)

[0074] Where Cov represents the covariance matrix, ave represents the mean temperature of transform_image(.) over all frames, and T represents the transpose. The eigenvectors of the covariance matrix Cov are the principal components with the largest variance in the PCA output. These principal components have only one dimension, which is the temperature value after removing the frame number.

[0075] The Cov eigenvectors of the covariance matrix are obtained using the `eig` function in MATLAB, denoted as `eig(Cov)`. These eigenvectors, representing the principal component components, are then multiplied by the physically reduced 2D data matrix `transform_image(.)` to obtain the PCA-transformed image data `PCAimage(·)`, which is the eigenvector of the Cov covariance matrix, representing the principal component with the largest variance in the PCA output. This allows for the next step of Independent Component Analysis (ICA).

[0076] PCAimage(h)=eig(Cov)transform_image(.) (3)

[0077] (2) Fixed-point iterative ICA and optimal solution;

[0078] For simplicity, let the image data after the PCA conversion be V:

[0079] V = PCAimage(h) (4)

[0080] (3) Update the weight matrix w by maximizing the negative entropy;

[0081] ICA improves data independence by enhancing the non-Gaussianity of the data, thus obtaining independent component units. To measure specific Gaussianity or non-Gaussianity, the concept of negative entropy is introduced; maximum negative entropy indicates maximum non-Gaussianity, meaning the data achieves the highest independence. W is calculated. T The negative differential entropy of V requires knowing W. T The probability density distribution function of V is complex to solve, so the approximate formula of equation (5) is usually used to conveniently solve for W. T The maximum negative differential entropy of V is shown in equation (5).

[0082] max(E{Vg(w T V)}-E{g′(w T V)}w)→w best (5)

[0083] In the formula, max(·) represents the maximum value of the polynomial; w best The weight matrix represents the value at which the negative entropy is maximized, indicating that the separation of independent components has been completed; g(·) is the nonlinear function for approximation, and g can be taken as y 3 ,tanh(a1y),yexp(-y 2 / 2); T represents matrix transpose, and E represents mathematical expectation.

[0084] Based on the optimal solution theory in nonlinear programming, the equality constraint optimization problem is solved through the Kuhn-Tucker condition, and the objective function for optimization is given by equation (5); under the constraint E{(w T V) 2}=||w|| 2 =1 (this constraint relationship is the inherent mathematical relationship after the data has been preprocessed), according to this constraint and objective function (5), after passing the Kuhn Tucker condition, we get equation (6), where β is the Kuhn Tucker non-negative coefficient.

[0085] E{Vg(w T V)}+βw=0 (6)

[0086] Newton's iteration method (Newton-Raphson method) is a commonly used method for solving complex equations. It is a method that sets an initial approximate solution and continuously approximates it through the intersection of tangents. The equation to be solved is equation (6), which is obtained as equation (7) through Newton's iteration method. Given the objective function (5) and the constraints of the Kuhn-Tucker condition, the weight matrix w is calculated using Newton's iteration method. k The optimization is performed using Newton's iteration method, as shown in equations (7) and (8), and the iteration stops when the ε threshold is met. k Set to w best .

[0087]

[0088] Among them, w k+1 This represents the iterative weight vector at step k+1.

[0089] (4) When the newly generated w k+1 Satisfy |w k+1 -w k |<ε, where ε is the threshold, usually set to 0.001, then let w best =w k+1 For this weight matrix w best Standardization processing, w best ←w best / ||w best ||;||.|| represents the 2-norm.

[0090] (5) Take the w obtained in step (4) bestMultiplying by V yields the independent image vector after whitening and ICA. The image vector processed by FASTICA can more prominently reflect the defect information of the detected object, that is, the image features are enhanced. Similarly, the dimension of the independent image vector is restored using the three-layer for loop program of MATLAB as described above, and it is restored from a one-dimensional vector to a two-dimensional vector (dimensionality increase). The two-dimensional vector is the number of rows and columns (i.e., the reverse processing of the above equation (1)). The specific two-dimensional vector value is temperature. The optimized image is as follows. Figure 5 As shown.

[0091] FASTICA derives the weight matrix within the specified error range using the Newton-Raphson iteration method, as shown in Equation (8). It then performs image processing using the opposite logic to dimensionality reduction (dimensionality increase), extracting key defect information by improving data independence. In data analysis, FASTICA can only consider improving data independence, requiring different nonlinear functions g(·) to process different images. For example, the cubic power function g works well for board 2, while board 1 requires the hyperbolic tangent function to achieve better image quality, with the hyperbolic tangent function coefficient set to 2. Furthermore, the output image data of the FASTICA processing section is set to I(x, y), where I(x, y) is the original image value.

[0092] I(x,y)=w best V (9)

[0093] Step 3: Multi-scale Retinex Image Denoising

[0094] In step 2, after FASTICA from equation (5) to equation (8), the maximum negative entropy weight w is obtained. best After multiplying with the image matrix V (i.e., equation (9)), the image is further subjected to multi-scale Retinex image denoising processing to obtain an infrared image that combines multi-scale Retinex and FASTICA. The detailed texture information of the image can usually be extracted by Retinex theory. Retinex theory is suitable for approximate quantitative processing of infrared thermal images. The workflow is as follows: Decompose an image I(x,y) into illumination component L(x,y) and image reflection component r(x,y) (this technique is well known to those skilled in the art). (x,y) is the image pixel position. In mathematical processing, the illumination component L(x,y) is usually decomposed into the temporal convolution of the original image value I(x,y) and the low-pass filter function (this technique is well known to those skilled in the art). The image reflection component r(x,y) can then be expanded into the following equation (10).

[0095] r(x,y)=lnI(x,y)-ln[I(x,y)*LPF(x,y)] (10)

[0096] In the formula, I(x,y) represents the original image value, and LPF(x,y) represents the low-pass filter function.

[0097] Low-pass filtering of images is usually implemented using a Gaussian kernel, as shown in equation (11).

[0098]

[0099] In equation (11), Gaussian(x,y) is the Gaussian function, λ is the normalization constant, and δ is the wrapping parameter of the Gaussian kernel. The larger the wrapping parameter δ, the larger the low-pass filtering range. The wrapping parameter δ is the scale of Retinex. Single scale means that the wrapping parameter δ has a unique value, and multi-scale means that the wrapping parameter δ has a non-unique value. By weighting r(x,y), based on the strong effectiveness of the Gaussian kernel's influence within 2 and 3 times the wrapping parameter δ, this invention sets three wrapping parameters for multi-scale Retinex processing: the first wrapping parameter δ1 is set to 128, the second wrapping parameter δ2 to 256, and the third wrapping parameter δ3 to 512. The low-pass filtering of multiple wrapping parameters in equation (12) is averaged (i.e., multi-scale). The image is filtered by low-pass filtering of different parameters, and the average value is taken to improve the low-pass filtering effect (i.e., using r). M (x, y) replaces equation (10) to improve the output illumination components. The multi-scale retinex expansion (12) is obtained from equation (11).

[0100]

[0101] Step 4: Signal-to-noise ratio calculation

[0102] To quantitatively compare the image optimization effects and verify the image processing method of this patent, this invention employs a comparative signal-to-noise ratio (SNR) calculation method. This method calculates the difference between the image vectors (the difference between two images) and the mean square error of all points in the image to obtain the SNR* value. When the SNR* value is greater than 1, the image quality is improved relative to the original image. Let the two image vectors being compared be Y and X, where the function num(Y) represents the number of elements in Y.

[0103] A = YX (13)

[0104] MSE=sum(A(:).*A(:)) / num(Y) (14)

[0105] SNR * =10*log 10 (sum(X(:).*X(:)) / MSE / num(Y)) (15)

[0106] Where A represents the multi-channel pixel difference map of the two images, sum(·) represents the summation of all points within the parentheses, A(:) represents all pixels in the pixel difference map, X(:) represents all pixel values ​​in the original image vector, .* represents matrix dot product, and SNR is calculated. * This indicates the signal-to-noise ratio.

[0107] After carbon fiber laminates No. 1 and No. 2 were uniformly heated by laser infrared, the FASTICA infrared thermographic processing method based on Retinex was applied to enhance the images of defects in the two carbon fiber laminates. Figure 2 , 5 The dimensions and defect distribution of plates 1 and 2 are shown respectively. Figure 3 , 4 Images 6 and 7 correspond to the processing effects of FASTICA and the combined Retinex image, respectively. Firstly, the corresponding defects are clearly visible, and these defects are more obvious compared to infrared thermography. The uniform light spot generated by laser heating is also well distinguishable. The two materials have the same number of Retinex scales and parameter values, but different FASTICA nonlinear functions.

[0108] When using the cubic nonlinear function, the SNR* value of board 1 was 6.9313. When the nonlinear function in FastICA was changed to the tanh function, the SNR* value was 14.5943. When board 2 changed to the cubic function, the SNR* value was 12.9048. It can be seen that the image processing effect was enhanced.

Claims

1. A multi-scale Retinex and Fastica combined infrared thermal imaging processing method, characterized in that, Includes the following steps: Step 1: Laser infrared thermal imaging; A laser infrared thermal imaging non-destructive testing system was used to perform laser infrared non-destructive testing on carbon fiber laminates, and a sequence of num frames of infrared thermal images was obtained through laser infrared thermal imaging. Step 2: Fastica image processing; (1) PCA processing of image data; The infrared thermal image sequence is whitened using the PCA preprocessing method, as shown in Equation (1). The infrared thermal image sequence acquired by the laser infrared thermal imaging non-destructive testing system is physically reduced in dimension, and the three-dimensional temperature data sequence is converted into a two-dimensional matrix composed of multiple one-dimensional data vectors. The three dimensions are the frame number, row number, and column number, respectively. The value of the element in the three-dimensional matrix determined by these three dimensions is the temperature value corresponding to the pixel at the frame number, row number, and column number. The three-dimensional temperature data sequence is converted into a two-dimensional matrix composed of multiple one-dimensional data vectors: transform_image(i,h)=image(i,j,v) (1) Where transform_image(.) represents a two-dimensional matrix composed of multiple one-dimensional data vectors of infrared thermograms. For two-dimensional data of the same frame, the second column of the original data is appended to the end of the first column, the third column is appended to the end of the second column, and so on. The two-dimensional data of the same frame becomes a one-dimensional data sequence. In equation (1), i represents the number of image frames, h represents the one-dimensional temperature sequence of a single frame after image dimensionality reduction, and image(.) represents a three-dimensional temperature data sequence. Where i, j, and v represent the number of frames, rows, and columns corresponding to a certain element in the three-dimensional temperature matrix, respectively. To facilitate data processing and fitting, mean reduction is required, which means removing the mean from the above two-dimensional matrix to obtain a zero-mean vector. After mean reduction and covariance calculation, the covariance matrix Cov is obtained. Cov=(transform_image(.)-ave)(transform_image(.)-ave) T (2) Where Cov represents the covariance matrix, ave represents the mean temperature of transform_image(.) over all frames, and T represents the transpose; the eigenvector of the covariance matrix Cov is the principal component with the largest variance in the PCA output. This principal component has only one dimension, which is the temperature value after removing the frame number. The Cov eigenvector of the covariance matrix is ​​obtained, denoted as eig(Cov). The eigenvector eig(Cov) representing the principal component is multiplied by the physical dimension-reduced two-dimensional data matrix transform_image(.) to obtain the PCA-transformed image data PCAimage(h), which is the eigenvector of the covariance matrix Cov mentioned above, that is, the principal component with the largest variance output by PCA. PCAimage(h)=eig(Cov)transform_image(.) (3) (2) Fixed-point iterative ICA and optimal solution; For simplicity, let the image data after the PCA conversion be V: V = PCAimage(h) (4) (3) Update the weight matrix w by maximizing the negative entropy. k ; Introducing the concept of negative entropy, the maximum negative entropy indicates the maximum non-Gaussianity, meaning the data achieves the highest independence; calculating W... T The negative differential entropy of V requires knowing W. T The probability density distribution function of V is complex to solve, so the approximate formula of equation (5) is used to solve for W. T The maximum negative differential entropy of V is shown in equation (5); max(E{Vg(w T V)}-E{g′(w T V)}w)→w best (5) In the formula, max(·) represents the maximum value of the polynomial; w best represents the weight matrix corresponding to the maximum negative entropy, indicating that the separation of independent components has been completed; g(·) is the nonlinear function for approximation; T represents matrix transpose, and E represents the expected value; Based on the optimal solution theory in nonlinear programming, the equality constraint optimization problem is solved through the Kuhn-Tucker condition, and the objective function for optimization is given by equation (5); under the constraint E{(w T V) 2 }=||w|| 2 =1, according to the constraints and objective function (5), the Kuhn Tucker condition is used to obtain equation (6), where β is the Kuhn Tucker non-negative coefficient; E{Vg(w T V)}+βw=0 (6) The equation to be solved is equation (6), which is then used to obtain equation (7) via Newton's iteration method. Given the objective function (5) and the constraints of the Kuhn-Tucker condition, the weight matrix w is calculated using Newton's iteration method. k The optimization is performed using Newton's iteration method, as shown in equations (7) and (8), and the iteration stops when the ε threshold is met. k Set to w best ; Among them, w k+1 This represents the iteration weight vector at step k+1; (4) When the newly generated w k+1 Satisfy |w k+1 -w k |<ε, where ε is the threshold, then let w best =w k+1 ; Regarding the weight matrix w best Standardization processing, w best ←w best / ||w best ||;||.|| denotes the 2-norm; (5) Take the w obtained in step (4) best Multiply by V to obtain the independent image vector after whitening and ICA; restore the dimension of the independent image vector, converting it from a one-dimensional vector to a two-dimensional vector, where the two dimensions are the number of rows and columns, and the specific two-dimensional vector value is temperature; FASTICA uses the Newton iteration method to obtain the weight matrix within the specified error, as shown in Equation (8). It uses the opposite program logic to dimensionality reduction, that is, the method of dimensionality increase, to process images. It extracts the key information of defects by improving data independence. In data analysis, FASTICA can only consider improving data independence. Different nonlinear functions g(·) are required to process different images. The output image data of the FASTICA processing part is set as I(x,y), where I(x,y) is the original image value. I(x,y)=w best V (9) Step 3: Multi-scale Retinex image denoising; In step 2, after FASTICA from equation (5) to equation (8), the maximum negative entropy weight w is obtained. best After multiplying with the image matrix V, the image is further subjected to multi-scale Retinex image denoising to obtain an infrared image combining multi-scale Retinex and FASTICA; the detailed texture information of the image is extracted by Retinex theory. The workflow is as follows: an image I(x,y) is decomposed into illumination component L(x,y) and image reflection component r(x,y), where (x,y) is the image pixel position. The illumination component L(x,y) is decomposed into the temporal convolution of the original image value I(x,y) and the low-pass filter function. The image reflection component r(x,y) can then be expanded into the following formula (10). r(x,y)=ln I(x,y)-ln[I(x,y)*LPF(x,y)] (10) In the formula, I(x,y) is the original image value, and LPF(x,y) represents the low-pass filter function; Low-pass filtering of images is usually implemented using a Gaussian kernel, as shown in equation (11); In equation (11), Gaussian(x,y) is the Gaussian function, λ is the normalization constant, and δ is the wrapping parameter of the Gaussian kernel. The larger the wrapping parameter δ, the larger the low-pass filtering range. The wrapping parameter δ is the scale of Retinex. Single scale means that the wrapping parameter δ has a unique value, and multi-scale means that the wrapping parameter δ has a non-unique value. By weighting r(x,y), based on the strong effectiveness of the Gaussian kernel's influence within 2-3 times the wrapping parameter δ, multiple wrapping parameters are set for multi-scale Retinex processing. The low-pass filtering of the image is performed using multiple wrapping parameters in equation (12). The average value of r(x,y) with different δ values ​​is taken to improve the low-pass filtering effect, that is, r is used to... M (x,y) replaces r(x,y) in equation (10), i.e., multi-scale, to improve the output illumination components; the multi-scale retinex expansion (12) is obtained by referring to equation (11). Where N represents the number of wrapping parameters, n represents the index of the wrapping parameter, and δ n Represents the nth wrap-around parameter, W n This represents the weighted value corresponding to the nth wrapping parameter; Step 4: Calculate the signal-to-noise ratio; A comparative signal-to-noise ratio (SNR) calculation method is employed, which calculates the mean square error of all points in the image by measuring the difference between image vectors, thereby obtaining the SNR. * Value, when the signal-to-noise ratio (SNR) * A value greater than 1 indicates an improvement in quality compared to the original image; let the two image vectors being compared be Y and X, where the function num(Y) represents the number of elements in Y; A = YX (13) MSE=sum(A(:).*A(:)) / num(Y) (14) SNR * =10*log 10 (sum(X(:).*X(:)) / MSE / num(Y)) (15) Where MSE represents mean squared error, A represents the multi-channel pixel difference map between the two images, sum(·) represents the summation of all points within the parentheses, A(:) represents all pixels in the pixel difference map, X(:) represents all pixel values ​​in the original image vector, .* represents matrix multiplication, and SNR is... * This indicates the signal-to-noise ratio.

2. The infrared thermal imaging processing method combining multi-scale Retinex and FASTICA as described in claim 1, characterized in that, g(·) is y 3 ,tanh(a1y),yexp(-y 2 One of the types of / 2).

3. The infrared thermal imaging processing method combining multi-scale Retinex and FASTICA as described in claim 1, characterized in that, ε is set to 0.

001.

4. The infrared thermal imaging processing method combining multi-scale Retinex and FASTICA as described in claim 1, characterized in that, In step 3, three wrapping parameters are set for the multi-scale Retinex processing.

5. The infrared thermal imaging processing method combining multi-scale Retinex and FASTICA as described in claim 4, characterized in that, In step 3, the first surround parameter δ1 is set to 128, the second surround parameter δ2 is set to 256, and the third surround parameter δ3 is set to 512.

Images