Generalized inverse matrix-exterior penalty function (GIM-EPF) algorithm for data processing of multi-wavelength pyrometer (MWP)

Jiafeng Liang; Li Dai; Sheng Chen; Weihong Gu; Bo Peng; Nannan Jiang; Wenlong Song; Jian Xing

doi:10.1364/OE.26.025706

1. Introduction

Multi-wavelength pyrometry (MWP) is useful for the measurement of high or ultra-high temperature targets whose temperatures cannot be measured by thermocouples [1–6]. It uses the intensity of the emitted thermal radiation in multiple spectral channels to directly retrieve true temperature and associated material spectral emissivity based on Planck's formula [7–13]. Compared with other non-contact temperature measurement methods, such as Holographic Interferometry [14], CARS [15], and TDLAS [16] method, the MWP method has the advantages of simplicity of instrumentation, minimal influence from external parameters, and direct inversion of Planck’s formula to obtain true surface temperature and surface spectral emissivity. It has been widely used in jet and rocket engine testing, plasma diagnostics, and other applications [17,18]. Data processing in MWP requires the solution of n equations based on Planck's law for n spectral channels; however, these n equations contain n + 1 unknowns, namely the n unknown spectral emissivities and the true temperature [19]. Before 2017, all the data processing methods for MWP assumed beforehand a model for the emissivity. The models including emissivity-wavelength models [20–23] and emissivity-temperature models [24–26]. Obviously, if the model is not consistent with reality, the accuracy of the parameters derived from the measurements will be greatly affected.

In 2017, in order to achieve true temperature inversion without any emissivity assumption, we proposed that the MWP data processing problem be converted into a constraint optimization problem. Gradient Projection (GP) and Internal Penalty Function (IPF) algorithms are standard approaches for solving this kind of problem [27]. The results demonstrate that the IPF algorithm is much better than GP in both accuracy and efficiency. IPF needs as input a reasonable emissivity range. However, the efficiency of the IPF approach is poor, and simulation results for 6 kinds of materials show that the average inversion time of IPF takes about 100 seconds, while the GP method requires more than 250 seconds. This makes it difficult to meet the requirement of real-time temperature measurement. In the present paper we develop and test a new Generalized Inverse Matrix-Exterior Penalty Function (GIM-EPF) data processing algorithm. The GIM-EPF algorithm greatly improves the processing speed for MWP data while ensuring good temperature accuracy. Additionally, the algorithm does not require any advanced knowledge of the emissivity.

2. Principles of algorithm

2.1 Reference temperature model

If MWP has n spectral channels, output signal V_i of the ith channel is defined as follows:

V_{i} = A_{λ_{i}} ε (λ_{i}, T) \frac{1}{λ_{i}^{5} (e^{C_{2} / λ_{i} T} - 1)} (i = 1, 2, ..., n)

where A_λi is a calibration factor which depends on wavelength. Detector sensitivity, geometry, absorption through windows, and the first radiation constant are all taken into account in the constant A_λi. ε(λ_i, T) is spectral emissivity of target at true temperature T; C₂ is second radiation constant; λ_i is effective wavelength of MWP at the ith channel; T is true temperature of target.

In the case of reference temperature T´ of blackbody, the ith channel output signal V_i´ is given by:

V_{i}^{'} = A_{λ_{i}} λ_{i}^{- 5} e^{- \frac{C_{2}}{λ_{i} T^{'}}}, [ε (λ_{i}, T) = 1.0 \begin{matrix} f o r & b l a c k b o d y \end{matrix}]

the ratio of Eq. (1) and Eq. (2) is described by:

\frac{V_{i}}{V_{i}^{'}} = ε (λ_{i}, T) e^{- \frac{C_{2}}{λ_{i} T}} e^{\frac{C_{2}}{λ_{i} T^{'}}}

This model only needs to measure output V_i´ of each channel under reference temperature T´, as long as reference temperature T´ is stable, it does not affect calculation results of target true temperature T and spectral emissivity ε(λ_i, T). The reference temperature model are simple and its condition is easy to meet, and error is less than other models such as brightness temperature model and calibration constant model, so it is widely used in data processing of MWP.

2.2 Combination of data processing and generalized inverse matrix of MWP

According to Eq. (3), after performing logarithmic operations on both sides of Eq. (3), so

\ln (\frac{V_{i}}{V_{i}^{'}}) - \frac{C_{2}}{λ_{i} T^{'}} = - \frac{C_{2}}{λ_{i} T} + \ln ε (λ_{i}, T)

Left of Eq. (4) is known, right is unknown. Order

y_{i} = \ln (\frac{V_{i}}{V_{i}^{'}}) - \frac{C_{2}}{λ_{i} T^{'}}

,

x_{i} = \ln ε (λ_{i}, T)

,

- \frac{C_{2}}{λ_{i} T} = a_{i} Z

,

a_{i} = - - \frac{C_{2}}{λ_{i}}

,

Z = \frac{1}{T}

.

So Eq. (4) is transformed as follows:

\begin{array}{l} y_{1} = x_{1} + a_{1} Z \\ y_{2} = x_{2} + a_{2} Z \\ ⋮ \\ y_{n} = x_{n} + a_{n} Z \end{array}

that is

[\begin{matrix} y_{1} \\ y_{2} \\ ⋮ \\ y_{n} \end{matrix}] = [\begin{matrix} 1 & 0 & \dots & 0 & a_{1} \\ 0 & 1 & 0 & 0 & a_{2} \\ ⋮ & 0 & 1 & ⋮ & ⋮ \\ 0 & 0 & \dots & 1 & a_{n} \end{matrix}] \times [\begin{matrix} x_{1} \\ x_{2} \\ ⋮ \\ x_{n} \\ Z \end{matrix}]

abbreviated is

Y = A X

where A is n × (n + 1) matrix, X is (n + 1) × n matrix. Equation (7) is an underdetermined equation, it has infinitely groups of solutions.

Because A is not a square matrix, there is no inverse matrix A⁻¹. In order to solve such problems, n × n matrix inverse is extended to n × m matrix inverse, named generalized inverse matrix. Because solution of Eq. (7) is not unique, we further restrict the generalized inverse matrix, that is, Moore-Penrose inverse matrix is used to guarantee the uniqueness of solution.

The matrix G that satisfies the following properties, which is defined as the Moore-Penrose inverse matrix of matrix A, and recorded A⁺, it is given as follows:

\begin{array}{l} A G A = A \\ G A G = G \\ {(A G)}^{H} = A G \\ {(G A)}^{H} = G A \end{array}

where A^H is a transposed conjugate matrix of A.

If matrix G satisfies all four formulas of Eq. (8) at the same time, then A⁺ is unique. For m × n matrix A, when rank r = n (that is, A is full column rank),

A^{+} = {(A^{H} A)}^{- 1} A^{H}

when rank r = m (that is, A is full row rank)

A^{+} = A^{H} {(A^{H} A)}^{- 1}

For data processing problem of MWP in Eq. (7), A is n × (n + 1) matrix, rank r = n (A is full row rank), so A⁺ is expressed by Eq. (10), and the unique solution for Eq. (7) is given:

X = A^{+} Y

Where solution X will meet

\min_{Y = A x} | | x | |

, min

‖ \cdot ‖

means Euclidean norm minimum. At this point, minimum distance between points means convergence. For emissivity, its range is limited from 0 to 1, so when ε_i is to

\min_{Y = A x} | | ε | |

, the solution is convergent and closer to the real value.

2.4 Generalized inverse matrix algorithm simulation for initial emissivity

According to generalized inverse matrix idea, simulation experiment is performed on six kinds of target materials labeled A-F, which are representative of different change trend of emissivity. The spectral emissivity of A-F targets is shown in Fig. 1 (blue ‘*’ market), true temperature is 1800K, reference temperature of blackbody is 1600K, 8 effective wavelengths are 0.4μm, 0.5μm, 0.6μm, 0.7μm, 0.8μm, 0.9μm, 1.0μm, 1.1μm, respectively.

Fig. 1 Emissivity simulation result based on directly generalized inverse matrix only.

Download Full Size | PDF

As can be seen from Fig. 1, trend of six target emissivity inversion results is basically the same as the actual distribution, but inversion values are all larger than true values, and the inversion values of partial emissivity are more than 1, which is not in accordance with the actual situation. It can be seen from Table 1 that the mean absolute error of temperature inversion results of six targets is 61.8K at 1800K, average relative error is 3.44%, and accuracy of temperature inversion results is not high. However, the emissivity inversion results based on generalized inverse matrix can be used as initial point of constrained optimization algorithm, so as to improve further accuracy and efficiency of inversion.

Table 1. Temperature Simulation by Generalized Inverse Matrix Only (True Temperature T = 1800K)

View Table | View all tables in this article

2.3 Combination of data processing and exterior penalty function of MWP

According to Eq. (3), ideally inversion temperature T of every spectral channel should be equal to true temperature if ε(λ_i, T) is known. Therefore, theoretical deviation of calculated temperature should tend to zero for every channel, so it is represented by the mathematical expression:

\sum_{i = 1}^{n} {[T_{i} - E (T_{i})]}^{2} = 0

where E(T_i) is the average value of inversion temperature of all channels. In the process of calculation, Eq. (12) can be converted into constrain optimization problem by minimized the deviation and made it infinitely close to zero, that is:

\min f (x) = \sum_{i = 1}^{n} {[T_{i} - E (T_{i})]}^{2} \to 0

After performing logarithmic operations on both sides of Eq. (3), it is given:

\ln (\frac{V_{i}}{V_{i}^{'}}) - \frac{C_{2}}{λ_{i} T^{'}} = - \frac{C_{2}}{λ_{i} T} + \ln ε (λ_{i}, T)

Order $x_{i} = \ln ε (λ_{i}, T)$ , $D_{i} = \frac{C_{2}}{λ_{i} T^{'}} - \ln (\frac{V_{i}}{V_{i}^{'}})$ ,

then inversion temperature T_i of every channel and E(T_i) can be expressed as follows:

T_{i} = \frac{C_{2}}{λ_{i}} \frac{1}{x_{i} + D_{i}}

E (T_{i}) = \frac{1}{n} \sum_{i = 1}^{n} \frac{C_{2}}{λ_{i}} \frac{1}{x_{i} + D_{i}}

minf(x) is derived from Eq. (15) and Eq. (16) according to Eq. (13), it is given:

\min f (x)= \sum_{i = 1}^{n} {[\frac{C_{2}}{λ_{i}} \frac{1}{x_{i} + D_{i}} - \frac{1}{n} \sum_{i = 1}^{n} \frac{C_{2}}{λ_{i}} \frac{1}{x_{i} + D_{i}}]}^{2}

Since emissivity range is [0, 1], that is 0≤ε(λ_i, T)≤1, the constraint condition is x_i<0.

The combination of Eq. (17) and x_i<0, then constraint optimization problem of MWP is defined as follows:

{\begin{matrix} \min f (x) = \sum_{i = 1}^{n} {[\frac{C_{2}}{λ_{i}} \frac{1}{x_{i} + D_{i}} - \frac{1}{n} \sum_{i = 1}^{n} \frac{C_{2}}{λ_{i}} \frac{1}{x_{i} + D_{i}}]}^{2} \\ x_{i} < 0 \end{matrix} .

This is a standard constraint optimization problem, Eq. (18) can be resolved by some constraint optimization algorithms such as Internal Penalty Function (IPF) method. However, from the results obtained from generalized inverse in Fig. 1, we can see that partial emissivity value obtained by Eq. (11) has been beyond the range of [0,1], IPF is not applicable. Therefore, exterior penalty function constraint algorithm must be used, emissivity starter value may be out of the [0, 1] range for it.

For inequality constrained optimization problem

{\begin{cases} \min f (x) \\ g_{i} (x) \geq 0, i = 1, 2, ..., k, k \in Z \end{cases}

penalty function is shown as followed

F (x) = f (x) + P \sum_{i = 1}^{k} m i n [(0, g_{i}^{2} (x))]

where P is penalty factor. Because Eq. (18) can be conformed into Eq. (19), it can be calculated using exterior penalty function method to process data of MWP.

Programming in MATLAB to achieve the exterior penalty function: minPF.

Call format: [x, minf] = minPF(f, x₀, g, c₁, p, var, eps)

Where x is independent variable; minf is the minimum value of objective function; f is objective function; x₀ is initial point; g is constraint function; c₁ is initial value of penalty factor; p is ratio coefficient of penalty factor; var is argument vector; eps is precision.

2.4 Parameters determination of minPF

Some parameters of minPF need to be determined before utilizing exterior penalty function algorithm.

(1) Initial point x₀
Initial point may be any point for exterior penalty function algorithm, but if the initial point is near the true value, calculation efficiency will be much improved. So calculation result from generalized inverse matrix is as the initial point according to Eq. (11), it would be an excellent initial emissivity, so blind choice of initial emissivity points can be avoided.
(2) Objective function f
From Eq. (15), we can see that penalty condition of exterior penalty function is g_j(x) ≥0, i = 1,2,…k, k ϵ Z, where g_i(x) is a constraint function, for Eq. (18), $x_{i} = \ln ε (λ_{i}, T) < 0$ is the constraint function. After many simulation experiment, emissivity is an idea range from 0.4 to 0.9, that is 0.4≤ε_i≤0.9, then −0.9≤x_i≤-0.1. Finally, inequality group ${\begin{cases} x_{i} \geq -0 .9 \\ - x_{i} \geq 0.1 \end{cases}$ , that is ${\begin{cases} x_{i} +0 .9 \geq 0 \\ - x_{i} - 0.1 \geq 0 \end{cases}$ is used to match constraint function of Eq. (15). Therefore, two matrices are given as follows:
$x_{i} = {[\begin{matrix} x_{1} & x_{2} & x_{3} & x_{4} & x_{5} & x_{6} & x_{7} & x_{8} \end{matrix}]}^{T}$ $b = {[\begin{matrix} 0.1 & - 0.9 & 0.1 & - 0.9 & 0.1 & - 0.9 & 0.1 & - 0.9 \end{matrix}]}^{T}$ $A = [\begin{matrix} - 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & - 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & - 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 1 & 0 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & - 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & - 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & - 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 1 & 0 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & - 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 1 & 0 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & - 1 \\ 0 & 0 & 0 & 0 & 0 & 0 & 0 & 1 \end{matrix}]$
Then $g_{i} (x) = A \times x_{i} - b$ is as constraint condition for data processing problem of MWP.
(3) Initial value c₁ and ratio coefficient p
Penalty factor P (Eq. (20) of exterior penalty function algorithm is controlled by initial value c₁ and ratio coefficient p, and it has effect on precision and efficiency of algorithm. To make the unconstrained problem close to the original constraint problem, we should choose the larger possible penalty factor, but in order to reduce the difficulty of solving the unconstrained problem, the smaller penalty factor should be selected, so how to select an appropriate penalty factor is very important. Generally, c₁ϵ[0.05, 0.5], p ϵ [2,5],in order to determine the optimal c₁ and p, following simulation data are analyzed for 6 targets in Fig. 1.
Firstly, keeping unchanged value for p (when p = 2), taking c₁ = 0.05, 0.1, 0.2, 0.3, 0.4, 0.5 separately to calculate temperature for A-F six targets.

After c₁ = 0.05, 0.1, 0.2, 0.3, 0.4, 0.5 was selected respectively, but c₁ has no effect on inversion results of temperature and emissivity, and average relative error is controlled at about 0.56% for temperature (Table 2), inversion trend of emissivity is also close to the actual value (Table 3), it showed a good inversion result from GIM-EPF algorithm.

Table 2. Relativity Error (%) of Temperature Simulation Results by GIM-EPF Algorithm to Different c₁, but Different c₁has the Same Results (p = 2, True Temperature T = 1800K)

View Table | View all tables in this article

Table 3. Emissivity Simulation Results by GIM-EPF to Different c₁ (p = 2, True Temperature T = 1800K)

View Table | View all tables in this article

Secondly, affection on running time from different c₁ is considered by GIM-EPF algorithm when p = 2.

Running time of GIM-EPF algorithm in the 6 targets is statistically analyzed. It is known from Table 4 that when c₁ = 0.05 is used, running time is shortest and average time of 6 targets is less than 5s.

Table 4. Running Time (s) by GIM-EPF Algorithm to Different c₁ (p = 2, True Temperature T = 1800K)

View Table | View all tables in this article

Thirdly, affection on minimum of target function minf = Eq. (17) from different c₁ is calculated by GIM-EPF algorithm when p = 2.

From Eq. (17), the closer mf is to zero, the closer the result is to optimal solution. Comparing value of A-F targets mf (in Table 5), mf is the nearest to zero when c₁ = 0.05.

Table 5. Minimum of Target Function minf = Eq. (17) to Different c₁ (p = 2, True Temperature T = 1800K)

View Table | View all tables in this article

Fourthly, keeping unchanged value of c₁ (when c₁ = 0.05), taking p = 2, 3, 4, 5 separately to calculate temperature and emissivity for A-F targets by GIM-EPF.

After p = 2, 3, 4, 5 was selected respectively, but p has no effect on inversion results of temperature and emissivity, and average relative error is controlled at about 0.56% for temperature (Table 6), inversion trend of emissivity is also close to the actual value (Table 7), it showed a good inversion result using GIM-EPF algorithm.

Table 6. Relativity Error (%) of Temperature Simulation Results by GIM-EPF Algorithm to Different p, but Different p has the Same Results (c₁ = 0.05, True Temperature T = 1800K)

View Table | View all tables in this article

Table 7. Emissivity Simulation Results by GIM-EPF to Different c₁ (p = 2, True Temperature T = 1800K)

View Table | View all tables in this article

Fifthly, running time from different p is considered by GIM-EPF algorithm when c₁ = 0.05.

Running time of GIM-EPF algorithm for 6 targets is statistically analyzed. It is known from Table 8 that running time is shortest and average time of 6 targets is less than 5s when p = 2 is used.

Table 8. Running Time (s) by GIM-EPF Algorithm to Different p (c₁ = 0.05, True Temperature T = 1800K)

View Table | View all tables in this article

Lastly, affection on minimum of target function minf = Eq. (17) from different p is calculated by GIM-EPF algorithm when c₁ = 0.05.

From Eq. (17) knowing that the closer mf is to zero, the closer the result is to the optimal solution. Comparing value of A-F targets mf (in Table 9), mf is the nearest to zero when p = 2.

Table 9. Minimum of Target Function minf = Eq. (17) to Different p (c₁ = 0.05, True Temperature T = 1800K)

View Table | View all tables in this article

To sum up, when c₁ = 0.05, p = 2, both emissivity and temperature inversion results can meet requirements very well with the shortest running time and the highest calculation precision, so initial value of penalty factor c₁ and ratio coefficient of penalty factor p of GIM-EFP algorithm parameters are c₁ = 0.05, p = 2.

3. Comparison between GIM-EPF algorithm and others

In order to verify accuracy, efficiency and noise resistance of GIM-EPF algorithm, GIM-EPF and other constraint optimization algorithms in [27], such as GP and IPF, are compared respectively under the condition of no noise and 5% random noise on V_i in Eq. (1).

3.1 Comparison results with no noise

From Fig. 2, it can be seen that accuracy of emissivity inversion of GIM-EPF algorithm is higher than that of GP algorithm, which is very closer to that of IPF algorithm. Even in the case of complex emissivity distribution of F material, it is superior to IPF algorithm, which shows that accuracy of emissivity inversion of GIM-EPF algorithm is high.

Fig. 2 Six targets emissivity simulation result based on different algorithms with no noise.

Download Full Size | PDF

It can be seen from Table 10 that temperature inversion accuracy of GIM-EPF algorithm is better than that of GP algorithm, which is equivalent to IPF algorithm, and average relative error is 0.56% for A-F materials at 1800K. Therefore, accurate inversion of temperature can be achieved based on GIM-EPF algorithm.

Table 10. Relativity Error (%) of Temperature Simulation Results with Different Algorithm (True Temperature T = 1800K)

View Table | View all tables in this article

As can be seen from Fig. 3, inversion efficiency is improved greatly by GIM-EPF algorithm, average running time is 4.0s (Operating environment: Intel(R) Core(TM)i5-4210U CPU@1.70GHz 2.40GHz), and average calculation time of IPF algorithm is 34.6s, that of GP is 561.5s, GIM-IPF efficiency is improved 8.7 times higher than that of IPF, and it is more than 140 times higher than that of GP algorithm. This is just because initial value obtained by using generalized inverse matrix idea is brought into EPF, which greatly reduces searching scope of emissivity, thus greatly improving computing efficiency. So if computer performance is better, on-line and high precision data processing of MWP can be realized by using GIM-EPF algorithm.

Fig. 3 Running time comparison to different algorithms with no noise for A-F materials.

Download Full Size | PDF

3.2 Comparison results with 5% random noise

It can be seen from Fig. 4 that in the case of 5% random noise on V_i in Eq. (1), GIM-EPF algorithm can still obtain a better inversion result of emissivity. The new algorithm has a much higher inversion accuracy than GP and IPF algorithm, especially for the E and F materials with large emissivity-wavelength fluctuation, GIM-EPF algorithm has the best results, it shown that GIM-EPF algorithm is good at anti-noise performance.

Fig. 4 Six targets emissivity simulation result based on different algorithms with 5% random noise on V_i.

Download Full Size | PDF

It can be seen from Table 11 that temperature inversion accuracy of GIM-EPF algorithm is better than that of GP and IPF algorithm with 5% random noise on V_i in Eq. (1), and average relative error is 0.56% for A-F materials at 1800K, especially for the E and F materials with large emissivity-wavelength fluctuation, GIM-EPF algorithm is superior to IPF, it indicated that new algorithm can be adapt to any material with good performance at anti-noise aspect. Therefore, novel GIM-EPF algorithm can be as the MWP data process method.

Table 11. Relativity Error (%) of Temperature Simulation Results with Different Algorithm (True Temperature T = 1800K, 5% Random Noise on V_i)

View Table | View all tables in this article

The most obvious advantage of GIM-EPF algorithm compared with other algorithms is that the calculation time is greatly reduced. From Fig. 5, it can be seen that in the case of 5% noise, average running time of GIM-EPF algorithm is only 4.34s, which is equivalent to the 4.01s without noise, but that of the IPF algorithm is 34.14s, and the GP algorithm is 241.3s. This proves that initial emissivity value based on generalized inverse matrix greatly improves searching efficiency of exterior penalty function algorithm, and it embodies excellent operational efficiency and good noise resistance ability for proposed GIM-EPF algorithm.

Fig. 5 Six materials running time based on different algorithms with 5% random noise on V_i

Download Full Size | PDF

4. Inversion of rocket nozzle temperature by MWP

In order to further verify practicability of proposed GIM-IPF algorithm, Prof. SUN measured surface temperature of a rocket nozzle using 8 wavelengths MWP, and gave spectral channel voltage values at reference temperature 2252K and 12 continuous measuring time points of the rocket nozzle in [24]. Using the proposed GIM-IPF algorithm, inversion temperature of each time point is calculated and compared with IPF algorithm.

The theoretical true temperature of the engine flame near the nozzle exit indicated by the rocket engine designer is 2490.0K. Figure 6 is result of temperature inversion of each measurement time node. It shows that GIM-EPF algorithm has a higher inversion accuracy. Compared with the IPF algorithm, although the new algorithm is slightly lower than IPF at some time points, the maximum relative error is only 0.65%, which reflects the higher temperature inversion accuracy than IPF. As can be seen from Fig. 7, GIM-EPF algorithm has a much shorter operation time than IPF algorithm. The average running time of GIM-EPF algorithm is only 3.4s, while the average running time of IPF algorithm reaches 210s. If hardware performance of a computer is improved, calculation time will be shorter. The new algorithm provides a theoretical basis for accurate on-line temperature measurement.

Fig. 6 Temperature comparison of two algorithms.

Download Full Size | PDF

Fig. 7 Running time comparison of two algorithms.

Download Full Size | PDF

5. Conclusions

A generalized inverse matrix and exterior penalty function (GIM-EPF) approach was integrated into the processing of multi–wavelength pyrometry data without requiring any advanced knowledge of the emissivity model and without the need to constrain the range of the emissivity beyond requiring that the ln ε ≤ 1. The proposed GIM-EPF algorithm eliminates the influence of emissivity on the inversion results to determine temperature. The simulation results from six materials show that the average relative error at 1800K temperature is 0.56%, a temperature inversion accuracy almost equal to IPF. A important advantage over the IPF is the inversion speed. The calculation speed of the GIM-EPF algorithm with good anti-noise performance is 8 times faster than that of IPF algorithm and is more than 100 times faster than that of Gradient Projection (GP) algorithm. Data from rocket plume temperature measurement confirms that the GIM-EPF algorithm has good inversion accuracy and efficiency. In summary, the GIM-EPF algorithm greatly improves the processing efficiency of MWP data while ensuring good temperature accuracy. The algorithm does not require any advanced knowledge of the emissivity. The GIM-EPF algorithm thus greatly increases the utility of MWP for measurement of high and ultra-high temperatures.

Funding

National Natural Science Foundation of China (NSFC) (61405045, 31470714); Fundamental Research Funds for the Central Universities (2572017DB04); “Double First Class” Talents Project Funding of NEFU 2016; Science and Technology Innovation Talents from Harbin City (RC2014QN009026).

References

1. I. Rancdarbord, G. Baudin, M. Genetier, D. Ramel, P. Vasseur, J. Legrand, and V. Pina, “Emission of gas and Al₂O₃ smoke in Gas-Al particle deflagration: experiments and emission modeling for explosive fireballs,” Int. J. Thermophys. 39(3), 152–160 (2018).

2. T. Matsumoto, T. Koizumi, Y. Kawakami, K. Okamoto, and M. Tomita, “Perfect blackbody radiation from a graphene nanostructure with application to high-temperature spectral emissivity measurements,” Opt. Express 21(25), 30964–30974 (2013). [CrossRef] [PubMed]

3. S. V. Onufriev, “Measuring the temperature of substances upon fast Heating with a current pulse,” Bull. Russ. Acad. Sci., Physics 82(4), 372–379 (2018). [CrossRef]

4. A. Bendada and M. Lamontagne, “A new infrared pyrometer for polymer temperature measurement during extrusion moulding,” Infrared Phys. Technol. 46(1–2), 11–15 (2004). [CrossRef]

5. F. Benedic, P. Bruno, and P. Pigeat, “Real-time optical monitoring of thin film growth by in situ pyrometry through multiple layers and effective media approximation modeling,” Appl. Phys. Lett. 90(13), 134104 (2007). [CrossRef]

6. A. Seifter and D. C. Swift, “Pyrometric measurement of the temperature of shocked molybdenum,” Phys. Rev. B 77(13), 761–768 (2007).

7. T. Fu, Z. Wang, and X. Cheng, “Temperature measurements of a diesel fuel combustion with multicolor pyrometry,” J. Heat Transfer 132(5), 051602 (2010). [CrossRef]

8. S. Liu, P. Farahmand, and R. Kovacevic, “Optical monitoring of high power direct diode laser cladding,” Opt. Laser Technol. 64(4), 363–376 (2014). [CrossRef]

9. T. Fu, H. Zhao, J. Zeng, M. Zhong, and C. Shi, “Two-color optical charge-coupled-device-based pyrometer using a two-peak filter,” Rev. Sci. Instrum. 81(12), 124903 (2010). [CrossRef] [PubMed]

10. D. Ng and G. Fralick, “Use of a multiwavelength pyrometer in several elevated temperature aerospace applications,” Rev. Sci. Instrum. 72(2), 1522–1530 (2001). [CrossRef]

11. V. Scharf and A. Katzir, “Four-band fiber-optic radiometry for determining the “true” temperature of gray bodies,” Appl. Phys. Lett. 77(19), 2955–2957 (2000). [CrossRef]

12. T. Fu, J. Liu, and A. Zong, “Radiation temperature measurement method for semitransparent materials using one-channel infrared pyrometer,” Appl. Opt. 53(29), 6830–6839 (2014). [CrossRef] [PubMed]

13. D. Svet and S. Sergeev, “Triple-wavelength pyrometer that measures true temperature,” Meas. Tech. 54(11), 1273–1275 (2012). [CrossRef]

14. J. S. Pérez-Huerta, T. Saucedo-Anaya, I. Moreno, D. Ariza-Flores, and B. Saucedo-Orozco, “Digital holographic interferometry applied to the investigation of ignition process,” Opt. Express 25(12), 13190–13198 (2017). [CrossRef] [PubMed]

15. D. Han, A. Satija, J. Kim, Y. Weng, J. Gore, and R. P. Lucht, “Dual-pump vibrational cars measurements of temperature and species concentrations in turbulent premixed flames with CO₂ addition,” Combust. Flame 181(3), 239–250 (2000).

16. G. Zhang, J. Liu, Z. Xu, Y. He, and R. Kan, “Characterization of temperature non-uniformity over a premixed CH 4 –air flame based on line-of-sight TDLAS,” Appl. Phys. B 122(1), 3–11 (2016). [CrossRef] [PubMed]

17. M. Liang, B. Sun, X. Sun, and J. Xie, “Development of a new fiber-optic multi-target multispectral pyrometer for achievable true temperature measurement of the solid rocket plume,” Measurement 95, 239–245 (2017). [CrossRef]

18. T. Fu, J. Liu, M. Duan, and S. Li, “Sub-pixel temperature measurements in hypersonic plasma jet environments using high speed multispectral pyrometry,” J. Heat Transfer 140(7), 1601–1609 (2018). [CrossRef]

19. P. Coates, “The least-square approach to multi-wavelength pyrometry,” High Temp. High Press. 20(4), 071601 (1988).

20. G. Gathers, “Analysis of multiwavelength pyrometry using nonlinear chi-square fits and Monte Carlo methods,” Int. J. Thermophys. 13(3), 539–554 (1992). [CrossRef]

21. M. A. Khan, C. Allemand, and T. W. Eagar, “Noncontact temperature measurement. I. Interpolation based techniques,” Rev. Sci. Instrum. 62(2), 392–402 (1991). [CrossRef]

22. H. Madura, M. Kastek, and T. Piatkowski, “Automatic compensation of emissivity in three-wavelength pyrometers,” Infrared Phys. Technol. 51(1), 1–8 (2007). [CrossRef]

23. P. Hagqvist, F. Sikström, A. Christiansson, and B. Lennartson, “Emissivity compensated spectral pyrometry for varying emissivity metallic measurands,” Meas. Sci. Technol. 25(2), 405–412 (2014).

24. X. G. Sun, G. B. Yuan, J. M. Dai, and Z. X. Chu, “Processing method of multi-wavelength pyrometer data for continuous temperature measurements,” Int. J. Thermophys. 26(4), 1255–1261 (2005). [CrossRef]

25. J. Xing, S. Cui, W. Qi, F. Zhang, X. Sun, and W. Sun, “A data processing algorithm for multi-wavelength pyrometry-which does not need to assume the emissivity model in advance,” Measurement 67(5), 92–98 (2015). [CrossRef]

26. J. Xing, R. S. Rana, and W. Gu, “Emissivity range constraints algorithm for multi-wavelength pyrometer (MWP),” Opt. Express 24(17), 19185–19194 (2016). [CrossRef] [PubMed]

27. J. Xing, B. Peng, Z. Ma, X. Guo, L. Dai, W. Gu, and W. Song, “Directly data processing algorithm for multi-wavelength pyrometer (MWP),” Opt. Express 25(24), 30560–30574 (2017). [CrossRef] [PubMed]

Targets	A	B	C	D	E	F
Inversion	1740.8	1723.4	1728.6	1752.2	1746.7	1737.3
Abs Error	59.2	76.6	71.4	47.8	53.3	62.7
Relative Error(%)	3.3	4.3	3.9	2.7	2.9	3.5

λ(μm)	A	B	C	D	E	F
0.4	0.79	0.51	0.48	0.69	0.74	0.48
0.5	0.71	0.53	0.58	0.63	0.58	0.54
0.6	0.64	0.55	0.68	0.57	0.50	0.64
0.7	0.58	0.59	0.78	0.49	0.59	0.82
0.8	0.53	0.64	0.77	0.49	0.78	0.64
0.9	0.48	0.71	0.67	0.59	0.61	0.54
1.0	0.47	0.79	0.57	0.69	0.52	0.64
1.1		0.44	0.46	0.74	0.47	0.84

c ₁	A	B	C	D	E	F	Average
0.05	4.228	4.914	3.479	6.459	5.601	4.976	4.943
0.1	4.212	6.443	6.848	4.742	3.526	5.257	5.171
0.2	3.978	7.379	6.708	4.915	4.43	6.286	5.616
0.3	6.474	5.616	5.554	6.615	6.615	6.895	6.295
0.4	4.914	4.789	7.784	6.224	6.459	4.181	5.725
0.5	5.335	6.583	6.115	5.148	6.38	6.349	5.985

c ₁	A	B	C	D	E	F	Average
0.05	2.46E-06	1.60E-06	2.12E-06	1.06E-06	1.60E-06	2.16E-06	1.83E-06
0.1	9.85E-06	6.39E-06	8.47E-06	4.22E-06	6.41E-06	8.64E-06	7.33E-06
0.2	3.94E-05	2.56E-05	3.39E-05	1.69E-05	2.56E-05	3.46E-05	2.93E-05
0.3	8.86E-05	5.75E-05	7.62E-05	3.80E-05	5.77E-05	7.77E-05	6.60E-05
0.4	1.57E-04	1.02E-04	1.35E-04	6.75E-05	1.03E-04	1.38E-04	1.17E-04
0.5	2.46E-04	1.60E-04	2.12E-04	1.06E-04	1.60E-04	2.16E-04	1.83E-04

λ(μm)	A	B	C	D	E	F
0.4	0.79	0.51	0.48	0.69	0.74	0.48
0.5	0.71	0.53	0.58	0.63	0.58	0.54
0.6	0.64	0.55	0.68	0.57	0.50	0.64
0.7	0.58	0.59	0.78	0.49	0.60	0.82
0.8	0.53	0.64	0.77	0.49	0.78	0.64
0.9	0.48	0.71	0.67	0.59	0.61	0.54
1.0	0.47	0.79	0.57	0.69	0.52	0.64
1.1	0.44	0.89	0.46	0.74	0.47	0.84

Generalized inverse matrix-exterior penalty function (GIM-EPF) algorithm for data processing of multi-wavelength pyrometer (MWP)

Abstract

1. Introduction

2. Principles of algorithm

2.1 Reference temperature model

2.2 Combination of data processing and generalized inverse matrix of MWP

2.4 Generalized inverse matrix algorithm simulation for initial emissivity

2.3 Combination of data processing and exterior penalty function of MWP

2.4 Parameters determination of minPF

3. Comparison between GIM-EPF algorithm and others

3.1 Comparison results with no noise

3.2 Comparison results with 5% random noise

4. Inversion of rocket nozzle temperature by MWP

5. Conclusions

Funding

References

Cited By

Figures (7)

Tables (11)

Equations (23)

Optics Express

Algorithm	A	B	C	D	E	F	Average
GP	0.81	3.84	3.84	0.81	0.81	3.35	2.24
IPF	0.26	0.05	0.38	0.95	0.78	0.42	0.47
GIM-EPF	0.37	0.67	0.36	1.10	0.72	0.16	0.56

Algorithm	A	B	C	D	E	F	Average
GP	1.61	0.60	1.18	2.49	2.39	1.47	1.62
IPF	0.72	0.22	0.41	0.93	0.79	0.44	0.59
GIM-EPF	0.35	0.66	0.42	1.11	0.75	0.09	0.56

p	A	B	C	D	E	F	Average
2	4.228	4.914	3.479	6.459	5.601	4.976	4.943
3	4.524	6.926	4.462	5.304	5.912	4.899	5.338
4	3.775	4.867	7.612	5.117	4.493	6.474	5.390
5	6.537	7.706	3.900	8.050	4.228	5.694	6.019

p	A	B	C	D	E	F	Average
2	4.228	4.914	3.479	6.459	5.601	4.976	4.943
3	4.524	6.926	4.462	5.304	5.912	4.899	5.338
4	3.775	4.867	7.612	5.117	4.493	6.474	5.390
5	6.537	7.706	3.900	8.050	4.228	5.694	6.019