Expand this Topic clickable element to expand a topic
Skip to content
Optica Publishing Group

Calibration-free wavelength modulation spectroscopy based on even-order harmonics

Open Access Open Access

Abstract

This paper proposes a novel and rapid calibration-free wavelength modulation spectroscopy algorithm based on even-order harmonics. The proposed algorithm, analytically deduced from Voigt line-shape function, only involves simple algebraic operations to describe the actual gas absorption spectra, thus eliminating the time-consuming simulations and line-shape fitting procedures adopted in traditional algorithms. Instead of acquiring the entirely scanned absorption line-shape, the proposed technique only requires extraction of the peak values of the harmonics. This characteristic significantly benefits gas diagnosis at elevated pressure and/or temperature, in which the entirely scanned absorption is very difficult to be obtained due to the broadened line-shapes. The proposed algorithm is validated by both numerical simulation and condition-controlled experiment, indicating millisecond-level calculation of gas parameters with the relative error less than 4% in the experiments.

© 2021 Optical Society of America under the terms of the OSA Open Access Publishing Agreement

1. Introduction

Tunable diode laser absorption spectroscopy (TDLAS) has been widely employed to measure gas properties, such as temperature and species concentration, given its advantages of contactless, rapid, and high sensitivity [15]. As a modality of the TDLAS technique, wavelength modulation spectroscopy (WMS) can achieve very high sensitivity with strong robustness to background noise, and therefore, has been extensively applied in harsh environments with, for example, strong turbulence, high pressure and/or temperature. Although WMS technique offers many benefits, quantitative WMS measurement is challenging since harmonic signals that are used for calculating the gas parameters strongly depend on the absorption line-shape. Consequently, the measurement results generally need to be calibrated with a known gas mixture, which is impractical in case of (a) poorly-known gas conditions in industrial applications, (b) difficulty in deploying a well-controlled calibration equipment, in particular, with high-temperature and high-pressure gases, and (c) optical inconsistency, e.g., transmission of the optical windows, between the calibration equipment and the target reactors. In the past decade, the above-mentioned challenges have stimulated the development of multiple calibration-free WMS methods. Implementation of the calibration-free WMS can be mainly categorized as model-dependent and model-independent methods.

The model-dependent method is to calculate the gas parameters using the analytical model of WMS, alongside the well-characterized laser parameters. Li et al. [6] pioneered the development of the model-dependent method by simulating WMS signals as a function of laser parameters and absorption spectra. To eliminate the influence of light intensity fluctuations, Rieker et al. [7] normalized the second-order harmonic signal of WMS by the first-order harmonic signal, i.e., WMS-2f/1f, and further included laser-specific tuning characteristics in the spectral-absorption model to infer gas properties. However, the retrieval accuracy of gas parameters using these methods highly depends on comprehensive characterization of both the line-shape parameters, e.g., collisional broadening coefficients, of the target transition, and the laser parameters, e.g., amplitudes and phases of laser intensity modulation and frequency modulation. In practice, these parameters measured in the laboratory can be very different from those measured on-site, due to the complex gas composition and wide ranges of temperature and pressure.

In contrast, the model-independent method, achieved by waveform fitting of the experimentally measured WMS-2f/1f (or WMS-nf/1f) signal with a simulated one, is more suitable for practical applications since less strict or even no requirement is placed on characterizing the line-shape parameters and laser parameters. This advantage significantly speeds up its application towards in situ diagnosis of reactive flows in harsh environments, for example, at elevated pressure and/or temperature [810], as the retrieval is independent of the prior information of line broadening. However, the waveform fitting requires time-consuming iterations until the residual is converged to an acceptable level. When the target absorption is interfered by nearby transitions, the result of the multi-parameter fitting cannot be always unique, as the iterations can be converged to the local optimum instead of the global one. To address these issues, Chen et al. [11] recently proposed to use the peak values of zero-, second-, and fourth-order harmonics, enabling the direct and fast measurement without complex simulations and curve fittings. However, Chen's algorithm is derived from the Lorentzian line-shape function, which is not applicable in low-pressure scenarios. Furthermore, the zero-order harmonic signal is easily polluted by low frequency noise, leading to reduced fidelity when the measurement is deployed at the presence of environmental noise owing to mechanical vibration, beam steering and window fouling.

Although higher order harmonics have lower signal strength [12,13], they can benefitted from larger signal-to-background ratios (SBRs) [1417] and stronger robustness at the extreme conditions where absorption spectra are significantly broadened and large modulation index is required. In this work, we propose, for the first time, a calibration-free WMS algorithm based on even-order harmonics (order ≥ 2). The amplitude of each harmonic signal, used for retrieving the gas properties, is analytically derived from the Voigt line-shape function, enabling not only speedy measurement down to milliseconds, but also a general suitability for various degrees of line-shape broadening. The rest of the paper is organized as follows: Section 2 details the fundamentals of the proposed algorithm. Numerical simulation and proof-of-concept experiment are carried out to validate the proposed algorithm in Section 3 and Section 4, respectively. The paper is finally summarized in Section 5.

2. Methodology

2.1 Theoretical background

In this subsection, some basics are presented to facilitate the derivation of the proposed method. WMS is implemented by imposing a high-frequency sinusoidal modulation with angular frequency ω on a slowly varying diode laser injection current. As a result, the instantaneous laser frequency v and the laser intensity It can be expressed as:

$$v = {v_0} + a\cos (\omega t),$$
$${I_t} = {\bar{I}_0}(1 + \sum\limits_{k = 1}^\infty {{i_k}\cos (k\omega t + {\psi _k})} ),$$
where v0 is the laser center frequency, a the modulation depth. $\bar{I}_0$ is the averaged laser intensity. ik is the k-th order Fourier coefficient of the laser intensity. ψk is the k-th order phase shift between intensity modulation and frequency modulation. For an isolated absorption line, the modulation index m is defined by 2a/λ, where λ represents the full width at half maximum (FWHM) of the absorption line.

According to the Beer-Lambert law, the spectroscopic absorbance can be written as:

$$\alpha (v) ={-} \ln ({I_t}/{I_0}) = A\varphi (v),$$
where It and I0 are the transmitted and incident laser intensities, respectively. A is the integral absorption area. φ(v) is the line-shape function, described by Voigt line-shape function, is a convolution of Lorentzian and Gaussian line-shape functions, noted as φL(v) and φG(v), respectively. φ(v) can be approximated by [18]
$$\varphi (v) = {c_L}{\varphi _L}(v) + {c_G}{\varphi _G}(v),$$
$${\varphi _L}(v) = 2/(\pi \lambda ){[1 + {((v - {v_0})/(\lambda /2))^2}]^{ - 1}},$$
$${\varphi _G}(v) = 2/\lambda \sqrt {\ln 2/\pi } \textrm{exp} [ - \ln 2{((v - {v_0})/(\lambda /2))^2}],$$
where cL and cG are the weights of the Lorentzian and Gaussian broadening coefficients, λL and λG, respectively. cL, cG and λ can be calculated by
$${c_L} = 0.6818817 + 0.6129331d - 0.1838439{d^2} - 0.1156844{d^3},$$
$${c_G} = 0.3246017 - 0.6182531d + 0.1768139{d^2} + 0.1210944{d^3},$$
$$\lambda = 0.5346{\lambda _L} + \sqrt {0.2166{\lambda _L}^2 + {\lambda _G}^2} ,$$
where the line-shape parameter d is defined as
$$d = ({\lambda _L} - {\lambda _G})/({\lambda _L} + {\lambda _G}).$$

A dimensionless parameter Δ = 2(v - v0)/λ is introduced to describe the relative offset of the laser center frequency, and therefore, φ(v) can be reformulated as

$$\varphi (v) = \frac{2}{{\pi \lambda }}[{c_L}L(\Delta ,m) + {c_G}\sqrt {\pi \ln 2} G(\Delta ,m)],$$
where L(Δ, m) and G(Δ, m) are the time-dependent peak normalized Lorentzian and Gaussian linear functions respectively, and are expressed as
$$L(\Delta ,m) = {[1 + {(\Delta + m\cos (\omega t))^2}]^{ - 1}},$$
$$G(\Delta ,m) = \textrm{exp} [ - \ln 2{(\Delta + m\cos (\omega t))^2}].$$

Following Arndt's model [14,19], the analytical expressions of the n-th harmonic of L(Δ, m) can be written as

$${L_n}(\Delta ) = \frac{{{\varepsilon _n}{i^n}}}{{2{m^n}}}\frac{{{{[\sqrt {{{(1 - i\Delta )}^2} + {m^2}} - (1 - i\Delta )]}^n}}}{{\sqrt {{{(1 - i\Delta )}^2} + {m^2}} }} + c.c.,$$
where c.c. denotes the conjugate complex number of the preceding term. ε0 = 1 and εn = 2 for n ≥ 1. Given no offset of the laser center frequency, i.e., Δ = 0, the amplitude of each harmonic obtained from Eq. (14) can be analytically expressed by
$$h_n^L = \left\{ {\begin{array}{cc} {\frac{{{\varepsilon_n}}}{{{m^n}}}\frac{{{{(\sqrt {1 + {m^2}} - 1)}^n}}}{{\sqrt {1 + {m^2}} }}}&{\textrm{for even }n}\\ 0&{\textrm{for odd }n} \end{array}} \right.. $$

By using infinite series theory [20], the amplitude of each harmonic of G(Δ, m) can be analytically expressed by

$$h_n^G = \left\{ {\begin{array}{cc} {2\textrm{exp} ( - {m^2}\ln 2/2){I_{{n / 2}}}({m^2}\ln 2/2)}&{\textrm{for even }n}\\ 0&{\textrm{for odd }n} \end{array}} \right.,$$
where In/2(z)=in/2Jn/2(iz). Jn/2 is the first kind of modified Bessel functions of order n/2.

As indicated by Eq. (11), φ(v) is a linear function of L(Δ, m) and G(Δ, m). When the linearity of Fourier transform is employed, the amplitude of the n-th harmonic hn of the given the spectral absorbance α(v) can be expressed by the following linear combination:

$${h_n} = \left\{ {\begin{array}{cc} {4A/\pi \lambda {c_L}{p_n}(m) + 4A/\pi \lambda {c_G}{q_n}(m)}&{\textrm{for even }n}\\ 0&{\textrm{for odd }n} \end{array}} \right.,$$
where pn(m) and qn(m) are defined as
$${p_n}(m) = \frac{{{\varepsilon _n}}}{{2{m^n}}}\frac{{{{(\sqrt {1 + {m^2}} - 1)}^n}}}{{\sqrt {1 + {m^2}} }},$$
and
$${q_n}(m) = \sqrt {\pi \ln 2} \textrm{exp} ( - {m^2}\ln 2/2){I_{{n / 2}}}({m^2}\ln 2/2),$$
respectively. Dependence of the amplitudes of pn(m) and qn(m) on the modulation index m is shown in Fig. 1. For a given m, both pn(m) and qn(m) decrease as the even harmonic order n increases. Consequently, the first several nonzero even harmonics should be adopted to achieve higher signal-to-noise ratios (SNRs) of the WMS measurement.

 figure: Fig. 1.

Fig. 1. Dependence of (a) pn(m) and (b) qn(m) on the modulation index m.

Download Full Size | PDF

2.2 Fundamentals of the proposed algorithm

The purpose of the proposed algorithm is to calculate the spectral parameters from the measured amplitudes of nonzero even harmonics. As indicated by Eq. (17), the amplitudes of the first N nonzero even harmonics adopted can be rewritten as:

$${\mathbf H} = {k_p}{\mathbf P}(m) + {k_q}{\mathbf Q}(m),$$
where H ∈ ℝN is the vector of measured harmonics amplitudes, defined as H= (h2, h4, …, h2N)T.

P(m) ∈ ℝN and Q(m) ∈ ℝN are m-dependent basis vectors, defined as P (m) = (p2(m), p4(m), …, p2N (m))T and Q(m) = (q2(m), q4(m), …, q2N(m))T, respectively. The coefficients kp and kq to be solved are defined as kp = 4AcL/(πλ) and kq = 4AcG/(πλ), respectively. Equation (20) shows that the vector H is an element of the column space in the matrix M = (P(m), Q(m)). Consequently, H and its projection on the column space of M, noted as H’, should coincide strictly:

$${\mathbf H} - {{\mathbf H}^\textrm{,}}{\mathbf = }({\mathbf E} - {\mathbf M}{({{\mathbf M}^\textrm{T}}{\mathbf M})^{ - 1}}{{\mathbf M}^\textrm{T}}){\mathbf H = 0},$$
where E is the N-dimensional identity matrix. Theoretically, with the measured H in hand, the modulation index can be calculated from the nonlinear Eq. (21). However, the measurement noise deviates the results of H - H’ in Eq. (21) from its desired value of 0. In this work, the measured modulation index m* is calculated by solving the following minimization problem:
$${m^\ast } = \mathop {\textrm{arg min}}\limits_m ||({\mathbf E - M}{({{\mathbf M}^\textrm{T}}{\mathbf M})^{ - 1}}{{\mathbf M}^\textrm{T}}){\mathbf H}|{|_2}.$$

When N = 2, M∈ ℝ2×2 is non-singular matrix, and thus the projection matrix, M(MTM)-1MT, is reduced to the identity matrix, E. In this case, H - H’ always equals to 0, and arbitrary modulation index m* can be chose. To make the above minimization problem solvable, the dimension of H must be no less than 3.

With the measured m* in hand, kp and kq in Eq. (20) can be calculated by using the least square method:

$$\left[ {\begin{array}{c} {{k_p}}\\ {{k_q}} \end{array}} \right]{\mathbf = }{({{\mathbf M}^{{\ast} \textrm{T}}}{{\mathbf M}^\ast })^{ - 1}}{{\mathbf M}^{\ast {\textrm T}}}{\mathbf H,}$$
where M* = (P(m*), Q(m*)). kp/kq is given by:
$$\frac{{{k_p}}}{{{k_q}}} = \frac{{{c_L}}}{{{c_G}}} = \frac{{0.6818817 + 0.6129331d - 0.1838439{d^2} - 0.1156844{d^3}}}{{0.3246017 - 0.6182531d + 0.1768139{d^2} + 0.1210944{d^3}}}.$$

Obviously, kp/kq is only determined by the line-shape parameter d. Figure 2 shows kp/kq is monotonically dependent on d when -0.85 ≤ d ≤ 0.85. That is to say, for any given d, there is a unique value of kp/kq. It is worth mentioning the line-shape is described as the Voigt function when -0.85 ≤ d ≤ 0.85. The line-shape functions can be characterized by the Lorentzian and the Gaussian functions for 0.85 < d ≤ 1 and -1 ≤ d < -0.85, respectively, resulting into simplified calculation of kp/kq to be detailed in subsection 2.3.

 figure: Fig. 2.

Fig. 2. Dependence of kp/kq on d when -0.85 ≤ d ≤ 0.85.

Download Full Size | PDF

Therefore, the FWHM of the Voigt line-shape can be calculated by

$$\lambda = 2a/{m^\ast },$$
and thus the integrated absorbance area can be calculated by
$$A = \pi \lambda {k_p}/(4{c_L})\textrm{ or }A = \pi \lambda {k_q}/(4{c_G}).$$

The integrated absorbance areas obtained from two individual transitions can be used to infer the gas temperature using two-line thermometry [21,22]. With the temperature in hand, the gas concentration can be calculated from one of the integrated absorbance areas. The flow chart for calculating gas properties using the proposed algorithm is shown in Fig. 3. It is worthy to point out that the logarithmic transformation method [23,24] was used to realize the normalization of harmonics in the proposed algorithm.

 figure: Fig. 3.

Fig. 3. Flow chart for calculating gas properties using the proposed algorithm.

Download Full Size | PDF

2.3 Simplified version for Lorentzian and Gaussian line-shape functions

As noted above, the line-shape function is Lorentzian when 0.85 < d ≤ 1, while that is Gaussian when -1 ≤ d < -0.85. In these two cases, Eq. (20) can be reformulated as

$${\mathbf H} = \left\{ {\begin{array}{cc} {{k_p}{\mathbf P}(m)}&{\textrm{for 0}\textrm{.85 < }d\textrm{ } \le \textrm{ }1\textrm{,}}\\ {{k_q}{\mathbf Q}(m)}&{\textrm{for - 1} \le d\textrm{ < - 0}\textrm{.85,}} \end{array}} \right.$$
where kp = kq = 4A/(πλ). Obviously, H is linear-dependent on P(m) or Q(m). Thus, the measured modulation index m* can be solved from the following minimization problems:
$${m^\ast } = \left\{ {\begin{array}{cc} {\mathop {\textrm{arg min }}\limits_m - \frac{{{\mathbf H} \cdot {\mathbf P}(m)}}{{||{\mathbf P}(m)|{|_2}}}}&{\textrm{for 0}\textrm{.85 < }d\textrm{ } \le \textrm{ }1\textrm{,}}\\ {\mathop {\textrm{arg min }}\limits_m - \frac{{{\mathbf H} \cdot {\mathbf Q}(m)}}{{||{\mathbf Q}(m)|{|_2}}}}&{\textrm{for - 1} \le d\textrm{ < - 0}\textrm{.85}\textrm{.}} \end{array}} \right.$$

Similarly, to make Eq. (28) solvable, the dimension of H must be no less than 2. Substituting m* into Eq. (27), kp and kq can be obtained as follows:

$${k_p} = \frac{{||{\mathbf H}|{|_2}}}{{||{\mathbf P}({m^\ast })|{|_2}}},$$
$${k_q} = \frac{{||{\mathbf H}|{|_2}}}{{||{\mathbf Q}({m^\ast })|{|_2}}}.$$

With kp or kq in hand, the line-shape parameters, such as FWHM, integrated absorbance area, and thus the gas parameters can be obtained in the same way described in subsection 2.2. Numerical results show that the deviation between the simplified version algorithm and the accurate proposed algorithm described in Section 2.2 is less than 0.113% when 0.85 < d ≤ 1 or -1 ≤ d < -0.85.

2.4 Parameter optimization

Equation (20) indicates H can be represented linearly by P(m) and Q(m). The coefficients kp and kq, therefore, are critical to solve the gas properties. According to the theory of linear algebra, the numerical stability for solving kp and kq can be characterized by the condition number of MTM. In case of an inverse problem with larger condition numbers, a small perturbation in H will cause large perturbations in kp and kq. Therefore, it is necessary to choose an appropriate modulation index m* to make kp and kq noise-insensitive.

Figure 4 shows the dependence of the condition number of MTM on modulation index m. It can be seen the condition number decreases as m starts to increase from 0, and reaches a local minimum value when m is around 1.5. Then, condition number increases and tops when m is around 3. After that, the condition number monotonically decreases as m increases. When m ≥ 3.5, the condition number is below the previous local minimum value, indicating the ever-stronger numerical stability for solving kp and kq. However, as indicated in Fig. 1, an excessively large m leads to the reduction of harmonic amplitude. Moreover, the wavelength tuning range of the laser is limited, which leads to the unreality of large modulation index. Therefore, the optimal m is selected to be 3.5.

 figure: Fig. 4.

Fig. 4. Dependence of the condition number of MTM on modulation index m.

Download Full Size | PDF

3. Numerical simulation

3.1 Simulation setup

The proposed algorithm is numerically validated for temperature measurement using water vapor (H2O) transitions. The absorption transitions at 7185.60 cm-1 and 7444.36 cm-1, that are widely used in two-line thermometry [21,22], are employed in both the simulation and experiment afterwards. In the simulation, the atmospheric-pressure gas concentration is set as 0.05. The absorption length is 20 cm. The gas temperature ranges from 500 K to 2,500 K. The frequencies of sinusoidal scan and modulation in the simulation are 100 Hz and 10 kHz, respectively. The modulation depths for the transitions at 7185.60 and 7444.36 cm-1 are set to 0.131 cm-1 and 0.127 cm-1, respectively. The scanning depths for the both transitions are set to 0.15 cm-1.

According to the HITRAN database [25], there are several neighbor transitions, with parameters detailed in Table 1, interfering the target transitions. Although each of the transitions appeared as a single absorption line-shape, the combined absorbance shown in Fig. 5 exhibits overlap of absorption line-shapes of nearby transitions with different lower energy levels. The overlap is particularly significant for the transition at 7444.36 cm-1, covering two interference transitions spaced up to 0.0133 cm-1. As shown in Fig. 5 (b), the two interference transitions result in significant asymmetry of the target absorbance with additional broadening of the line shape.

 figure: Fig. 5.

Fig. 5. The normalized line-shapes at different temperatures for the transitions at around (a) 7185.60 cm-1 and (b) 7444.36 cm-1, respectively.

Download Full Size | PDF

Tables Icon

Table 1. Parameters of the selected transitions at around 7185.60 cm-1 and 7444.36 cm-1. For each transition, the table shows its wavenumbers, line strengths (S(T0)) at T0 = 296 K, pressure broadening-induced half-widths for itself (ξself) and air (ξair), energies of the lower-state (E′′), and coefficients of temperature dependent broadening (nair) [25].

With the above-mentioned settings of the WMS simulation, the modulated absorbance signal at the temperature of 1000 K for the two target transitions are obtained, and thus the corresponding second- to tenth-order even harmonics are acquired by the FFT spectral analysis [26]. As shown in Fig. 6, the higher-harmonics give lower amplitude, but have the same order of intensities with the lower-order ones. In the proposed algorithm, only the peak values of the even harmonics are needed, as illustrated in Eqs. (20)–(30), thus significantly simplifying and fastening of data processing in comparison with the traditional waveform fitting algorithm.

 figure: Fig. 6.

Fig. 6. Simulated WMS raw signal and the extracted even harmonics acquired by the FFT spectral analysis at 1000 K. (a) and (b) show the modulated absorbance signal for the transitions at around 7185.60 cm-1 and 7444.36 cm-1 with the simulated H2O spectral absorption, respectively. (c) and (d) shows the second- to tenth-order even harmonics extracted from (a) and (b), respectively.

Download Full Size | PDF

3.2 Calculation of the modulation index

As noted in Eq. (22), the modulation index m* should be firstly solved by minimizing the objective function F = ||(E-M(MTM)-1MT)H||2, which depends on the selected orders of the even harmonics. For the two H2O transitions at 7185.60 cm-1 and 7444.36 cm-1, Fig. 7 shows the dependence of F on the modulation index m in three cases with different combinations of the nth order even harmonics, i.e., n = 2, 4, 6 (case 1), n = 2, 4, 6, 8 (case 2), and n = 2, 4, 6, 8, 10 (case 3), respectively. For both the transitions, F vs. m variation tendencies are same for all the three cases, with the minimum F obtained with similar m values. For the transitions at 7185.60 cm-1, the minimum F is obtained when m equals to 3.843, 3.866, and 3.872 for cases 1, 2 and 3, respectively. For the transitions at 7444.36 cm-1, the corresponding m equals to 3.796, 3.814, and 3.815 for cases 1, 2 and 3, respectively.

 figure: Fig. 7.

Fig. 7. Dependence of the objective function F on the modulation index m at 1000 K for the transitions at (a) 7185.60 cm-1 and (b) 7444.36 cm-1.

Download Full Size | PDF

3.3 Simulation results with noise-free and noise-contaminated data

With the modulation index m* and the amplitudes of the even-order harmonics in hand, we can calculate the integrated absorbance areas for the two transitions according to Eq. (26), and their ratio to infer the gas temperature. We will evaluate retrieved temperature using the proposed algorithm with noise-free and noise-contaminated measurements, respectively.

For the above three cases, Fig. 8 shows the relative error of the retrieved temperature using noise-free data. As the target transitions are interfered by the nearby transitions, minor errors exist in the retrieved temperature even if no noise is imposed on the measurement. For all the three cases, the relative errors are less than 0.5% in the temperature range of 500 - 1250 K. When the temperature is above 1250 K, the retrievals degenerate in all the three cases, with the errors in cases 1 and 2 less significant than that in case 3. The worst scenarios at 2500 K give the relative errors of -2.135%, -2.406%, and -3.175% for cases 1, 2, and 3, respectively.

 figure: Fig. 8.

Fig. 8. Relative errors of the retrieved temperature with different harmonic combinations.

Download Full Size | PDF

Then, white noise is imposed on both the transmitted and incident laser intensities, i.e., It and I0, according to the following formula:

$$\alpha ^{\prime}(v) ={-} \ln (\frac{{{I_0}\textrm{exp} ( - \alpha (v)) + nois{e_t}}}{{{I_0} + nois{e_0}}}),$$
where noiset and noise0 are white noise generated by the randn function of MATLAB with a mean value of 0 and a standard deviation of 0.01% × average laser intensity, equivalent to a signal-to-noise ratio (SNR) of 80 dB. This noise level is dominated by the thermal noise from the detector, signal amplification, and digitization circuits [27]. To characterize the noise performance of the proposed algorithm, the average relative error and the relative standard deviation (RSD) of the temperature retrieved from 100 repetitive simulations are shown in Figs. 9 (a) and (b), respectively. The average relative error shown in Fig. 9 (a) indicates similar results as those using the noise-free data. As the temperature sensitivities of the chosen H2O transitions decrease as the temperature increases [22,28], the retrieved temperature suffers from more significant fluctuations, indicated by the larger values of RSDs shown in Fig. 9 (b), at higher temperatures. In addition, involvement of more harmonics contributes to lower RSDs, and thus stronger robustness of the proposed algorithm.

 figure: Fig. 9.

Fig. 9. Performance evaluation of the proposed algorithm with noise-contaminated data. (a) and (b) show the average relative error and relative standard deviation (RSD) of the temperature retrieved from 100 repetitive simulations.

Download Full Size | PDF

4. Experimental setup and results

4.1 Experimental setup

In this section, experiments, with setup shown in Fig. 10, were carried out to validate the proposed algorithm. The quartz optical cell has three zones, all fitted with wedged windows to suppress etalon interference. The two outer zones were purged with dry air to prevent unwanted H2O absorption along the optical path, while the center zone, with the path length of 20 cm, was filled with a uniform gas mixture at a controlled temperature between 773 and 1273 K. The gas mixture is obtained from the dry air going through H2O liquid, and therefore, contains the target H2O molecules for the absorption measurement. Three K-type thermocouples were equally spaced along the optical cell to measure the temperature of the target gas. The optics, detector and the intermediate open paths were purged with N2 to remove interfering H2O absorption in the room air.

 figure: Fig. 10.

Fig. 10. Schematic diagram of the experimental system.

Download Full Size | PDF

Two distributed feedback (DFB) diode lasers operating near 7185.60 cm-1 (NTT, NLK1E5GAAA) and 7444.36 cm-1 (NTT, NLK1B5EAAA) were used to probe the H2O transitions. Each laser was controlled independently by a laser controller (Stanford Research Systems, LDC501), which adjusts the temperature and the injection current of the laser. To achieve wavelength modulation, the injection current was driven by a function generator (RIGOL, DG1062Z) which provided a 100 Hz sinusoidal scan signal superimposed with a 10 kHz sinusoidal modulation signal. The two lasers were operated in time division multiplexing (TDM) mode [29]. To monitor DC background signals such as the dark current of the detector and thermal radiation, there is a 10 ms non-laser interval between the neighboring wavelength scan of the two lasers, consequently, the overall temporal resolution of the measurements is 50 ms. The two lasers were combined and then split by a 2×2 single-mode fiber-coupler, with one output being guided through an etalon with a free spectral range (FSR) of 0.01 cm−1 for wavelength monitoring and the other output penetrating the target gas cell. Each laser beam was detected by a photodetector (Thorlabs, PDA10CS-EC) and further sampled by a 14-bit data acquisition card (NI Corporation, PXIe-5170R). The digitization rates for sampling the etalon signal and the absorption signal are 25 Mega samples per second (Msps) and 1 Msps, respectively.

The optical cell was firstly heated to the desired temperature and filled with dry air to obtain the absorption-free transmission. Then, the center zone was evacuated and filled with the H2O/air mixture. When the adsorption and transient effects was stabilized, the absorption transmissions were recorded at the given furnace temperature. The temperature of the optical cell was then adjusted to repeat the entire measurement process. In the experiment, the flow rates for both the purging gas and the target gas mixture were small and stabilized to reduce the influence of convective heat transfer on the temperature set by the furnace. Therefore, two mass flow meters were used to control both gas flow rates at 100 standard cubic centimeter per minute (sccm). Since the optical cell is open to the air, the gas pressure inside the cell is same as the atmospheric pressure, that is 0.998 atm, measured by the pressure gauge (MKS, type660).

4.2 Experimental results and discussions

Figures 11 (a) and (b) shows the sampled background transmission I0(t) and the absorption transmission It(t) at 973 K, respectively. Both I0(t) and It(t) are averaged for 100 repetitive measurements to reduce the interference of measurement noise. Although the background transmission is not an ideal absorption-free signal, caused by the small purging flow of the dry air, logarithmic transformation in Eq. (3) neutralizes any impact of the background absorption on the temperature retrieval. In addition, the rising edges of both the second periods of I0(t) and It(t) are adopted in the harmonic analysis in the experiment, since the first period, just after switching the lasers, is under stabilization. Subsequently, the even-order harmonics demodulated by FFT spectral analysis are shown in Figs. 12 (a) and (b) for the transitions at 7185.60 cm-1 and 7444.36 cm-1, respectively. Finally, the amplitude of each harmonic is extracted by conic fitting of the partial waveform nearby the peak.

 figure: Fig. 11.

Fig. 11. The experimentally sampled (a) background intensity and (b) transmitted intensity signals at 973 K.

Download Full Size | PDF

 figure: Fig. 12.

Fig. 12. The even-order (2nd to 10th) harmonics demodulated by FFT spectral analysis for the transitions at (a) 7185.60 cm-1 and (b) 7444.36 cm-1.

Download Full Size | PDF

With the optimal modulation index determined in subsection 3.2 and the extracted amplitudes of the even-order harmonics, the temperature is retrieved using the proposed algorithm by following Eqs. (20)–(30). To examine the performance of the proposed algorithm, the retrieved temperature values are compared with those measured by the thermocouple in the range of 773 - 1273 K with an interval of 100 K. Similar as the simulation, the three cases with different combinations of the even-order harmonics are verified in the experiment. As shown in Fig. 13(a), the temperature retrieved in all the three cases agree well with that measured by the thermocouple. The maximum relative errors, characterized in Fig. 13(b), are 2.59%, 3.60% and 2.82% for cases 1, 2 and 3, respectively.

 figure: Fig. 13.

Fig. 13. Comparison of the temperature measured using the proposed algorithm and the thermocouple. (a) and (b) show the absolute measured temperature and the relative error, respectively.

Download Full Size | PDF

To evaluate the superiority of the proposed algorithm in terms of computational cost, the waveform fitting algorithm [30] was also employed to analyze the same acquired data. Figure 14 shows the measured and fitted 1f-normalized WMS-2f (WMS-2f/1f) waveforms with the fitting residuals for the transitions of 7185.60 cm-1 and 7444.36 cm-1 at 973 K. The waveform fitting is implemented using the samples on the wavenumber free from the interference transitions. Therefore, 4000 samples on the left wing of absorbance are used for fitting for the transition 7185.60 cm-1 (Fig. 14(a)), while 2800 samples of the right wing are used for fitting for the transition 7444.36 cm−1.

 figure: Fig. 14.

Fig. 14. Measured and fitted WMS-2f/1f waveform for the transitions at (a) 7185.60 cm-1 and (b) 7444.36 cm-1 with the residuals shown in above subplots.

Download Full Size | PDF

It is worth noting the WMS-2f/1f waveform fitting algorithm is sensitive to the initial value of parameters to be fitted. Improper initial values may increase the computational time, and even make the solution fall into the local optimum. In this work, the furnace temperature given by thermocouples and other needed line-shape parameters retrieved by the proposed algorithm were taken as the initial values of the fitting algorithm. In contrast, the proposed algorithm is always converged with no restriction attached to the initial value. All the computational process in this section were conducted by a computer with an Intel Core i7-1165G7 @ 2.80 GHz CPU.

The computational cost results for these two algorithms are summarized in Table 2. The proposed algorithm in the temperature range of 773 - 1273 K gives the average computational time of 0.137 s, 0.142 s, and 0.144 s for cases 1, 2, and 3, respectively. Although the appropriate initial values were given, the average computational time of the WMS-2f/1f waveform fitting algorithm still reaches the average time of 22.9 s. That is to say, the computational efficiency of the proposed algorithm is improved by at least two orders of magnitude (159 - 167 times speed up).

Tables Icon

Table 2. Summary of computational time with different algorithms.

5. Conclusions

A novel calibration-free WMS algorithm, analytically deduced from the Voigt line-shape function, is proposed for accurate retrieval of gas properties using the amplitudes of the even-order harmonics. As only the peak values of the harmonics are needed, the proposed algorithm substantially benefits WMS measurement with large line-shape broadening, for example, measurements in high-pressure and high-temperature shock-wave wind tunnels. For the two pre-selected H2O transitions, numerical simulation was firstly carried out to determine the optimal modulation index for different combinations of the even-order harmonics, and examine their accuracy and noise resistance. Simulation results indicate the relative errors are less than 0.5% in the temperature range of 500 - 1250 K and stronger robustness when more harmonics are involved. Furthermore, experimental validation of proposed algorithm shows a maximum relative error of 3.6% in the temperature range of 773 - 1273 K. Compared with the traditional WMS-2f/1f waveform fitting algorithm, the computational efficiency achieved by the proposed algorithm is improved by at least two orders of magnitude. Finally, we would like to point out that the proposed algorithm is built upon the approximate description of the Voigt function given by the Eqs. (4)–(10). Therefore, the accuracy of the retrieved parameters cannot be better than the accuracy of this approximation. The updated algorithm built upon much more accurate Voigt approximation would be one of our future works.

Funding

National Key Research and Development Program of China (2017YFB0603204); Global Challenges Research Fund (LMIC Travel and Partnerships Fund).

Acknowledgement

Our deepest gratitude goes to the anonymous reviewers for their careful work and thoughtful suggestions that have helped improve this paper substantially.

Disclosures

The authors declare no conflicts of interest.

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

References

1. R. K. Hanson, R. M. Spearrin, and C. S. Goldenstein, Spectroscopy and optical diagnostics for gases (Springer, 2016).

2. R. K. Hanson, “Applications of quantitative laser sensors to kinetics, propulsion and practical energy systems,” Proc. Combust. Inst. 33(1), 1–40 (2011). [CrossRef]  

3. C. Liu and L. Xu, “Laser absorption spectroscopy for combustion diagnosis in reactive flows: A review,” Appl. Spectrosc. Rev. 54(1), 1–44 (2019). [CrossRef]  

4. W. Cai and C. F. Kaminski, “Tomographic absorption spectroscopy for the study of gas dynamics and reactive flows,” Prog. Energy Combust. Sci. 59, 1–31 (2017). [CrossRef]  

5. Y. Du, Z. Peng, and Y. Ding, “A high-accurate and universal method to characterize the relative wavelength response (RWR) in wavelength modulation spectroscopy (WMS),” Opt. Express 28(3), 3482–3494 (2020). [CrossRef]  

6. H. Li, G. B. Rieker, X. Liu, J. B. Jeffries, and R. K. Hanson, “Extension of wavelength-modulation spectroscopy to large modulation depth for diode laser absorption measurements in high-pressure gases,” Appl. Opt. 45(5), 1052–1061 (2006). [CrossRef]  

7. G. B. Rieker, J. B. Jeffries, and R. K. Hanson, “Calibration-free wavelength-modulation spectroscopy for measurements of gas temperature and concentration in harsh environments,” Appl. Opt. 48(29), 5546–5560 (2009). [CrossRef]  

8. I. A. Schultz, C. S. Goldenstein, C. L. Strand, J. B. Jeffries, R. K. Hanson, and C. P. Goyne, “Hypersonic scramjet testing via diode laser absorption in a reflected shock tunnel,” Journal of Propulsion Power 30(6), 1586–1594 (2014). [CrossRef]  

9. C. S. Goldenstein, R. M. Spearrin, and R. K. Hanson, “Fiber-coupled diode-laser sensors for calibration-free stand-off measurements of gas temperature, pressure, and composition,” Appl. Opt. 55(3), 479–484 (2016). [CrossRef]  

10. C. Liu, L. Xu, F. Li, Z. Cao, S. A. Tsekenis, and H. McCann, “Resolution-doubled one-dimensional wavelength modulation spectroscopy tomography for flame flatness validation of a flat-flame burner,” Appl. Phys. B 120(3), 407–416 (2015). [CrossRef]  

11. W. Chen, Y. Wu, and Y. Luo, “Broadening-independent demodulation for wavelength modulation spectroscopy based on even-order harmonics,” Spectrosc. Lett. 51(1), 61–66 (2018). [CrossRef]  

12. A. Upadhyay, D. Wilson, M. Lengden, A. L. Chakraborty, G. Stewart, and W. Johnstone, “Calibration-free WMS using a cw-DFB-QCL, a VCSEL, and an edge-emitting DFB laser with in-situ real-time laser parameter characterization,” IEEE Photonics J. 9(2), 1–17 (2017). [CrossRef]  

13. K. Sun, X. Chao, R. Sur, J. Jeffries, and R. Hanson, “Wavelength modulation diode laser absorption spectroscopy for high-pressure gas sensing,” Appl. Phys. B 110(4), 497–508 (2013). [CrossRef]  

14. S. Schilt, L. Thevenaz, and P. Robert, “Wavelength modulation spectroscopy: combined frequency and intensity laser modulation,” Appl. Opt. 42(33), 6728–6738 (2003). [CrossRef]  

15. J. Gustafsson, N. Chekalin, and O. Axner, “Characterization of 2f-, 4f-, and 6f-background signals in wavelength modulation diode laser absorption spectrometry in graphite furnaces,” Spectrochimica Acta Part B: Atomic Spectroscopy 58(1), 123–141 (2003). [CrossRef]  

16. P. Kluczynski and O. Axner, “Theoretical description based on Fourier analysis of wavelength-modulation spectrometry in terms of analytical and background signals,” Appl. Opt. 38(27), 5803–5815 (1999). [CrossRef]  

17. P. Kluczynski, Å. M. Lindberg, and O. Axner, “Background signals in wavelength-modulation spectrometry with frequency-doubled diode-laser light. I. Theory,” Appl. Opt. 40(6), 783–793 (2001). [CrossRef]  

18. Y. Liu, J. Lin, G. Huang, Y. Guo, and C. Duan, “Simple empirical analytical approximation to the Voigt profile,” J. Opt. Soc. Am. B 18(5), 666–672 (2001). [CrossRef]  

19. R. Arndt, “Analytical line shapes for Lorentzian signals broadened by modulation,” J. Appl. Phys. 36(8), 2522–2524 (1965). [CrossRef]  

20. P. Kluczynski, Å. M. Lindberg, and O. Axner, “Wavelength modulation diode laser absorption signals from Doppler broadened absorption profiles,” Journal of Quantitative Spectroscopy and Radiative Transfer 83(3-4), 345–360 (2004). [CrossRef]  

21. C. Liu, L. Xu, J. Chen, Z. Cao, Y. Lin, and W. Cai, “Development of a fan-beam TDLAS-based tomographic sensor for rapid imaging of temperature and gas concentration,” Opt. Express 23(17), 22494–22511 (2015). [CrossRef]  

22. C. Liu, L. Xu, Z. Cao, and H. McCann, “Reconstruction of Axisymmetric Temperature and Gas Concentration Distributions by Combining Fan-Beam TDLAS With Onion-Peeling Deconvolution,” IEEE Trans. Instrum. Meas. 63(12), 3067–3075 (2014). [CrossRef]  

23. Y. Wang, H. Cai, J. Geng, and Z. Fang, “Logarithmic conversion of absorption detection in wavelength modulation spectroscopy with a current-modulated diode laser,” Appl. Opt. 48(21), 4068–4076 (2009). [CrossRef]  

24. M. Cong and D. Sun, “Detection of harmonics and recovery of the absorption line profile using logarithmic-transformed wavelength modulation spectroscopy,” Opt. Eng. 55(7), 076119 (2016). [CrossRef]  

25. I.E. Gordon, L.S. Rothman, C. Hill, R.V. Kochanov, Y. Tan, P.F. Bernath, M. Birk, V. Boudon, A. Campargue, K.V. Chance, B.J. Drouin, J.-M. Flaud, R.R. Gamache, J.T. Hodges, D. Jacquemart, V.I. Perevalov, A. Perrin, K.P. Shine, M.-A.H. Smith, J. Tennyson, G.C. Toon, H. Tran, V.G. Tyuterev, A. Barbe, A.G. Császár, V.M. Devi, T. Furtenbacher, J.J. Harrison, J.-M. Hartmann, A. Jolly, T.J. Johnson, T. Karman, I. Kleiner, A.A. Kyuberis, J. Loos, O.M. Lyulin, S.T. Massie, S.N. Mikhailenko, N. Moazzen-Ahmadi, H.S.P. Müller, O.V. Naumenko, A.V. Nikitin, O.L. Polyansky, M. R. M. Rey, S.W. Sharpe, E. S. K. Sung, S.A. Tashkun, J. Vander Auwera, G. Wagner, J. Wilzewski, P. Wcisło, S. Yu, and E.J. Zak, “The HITRAN2017 molecular spectroscopic database,” Journal of Quantitative Spectroscopy and Radiative Transfer 203, 3–69 (2017). [CrossRef]  

26. S. Yin, Y. Liu, and T. Dong, “Trace gas detect based on spectral analysis and harmonic ratio,” Opt. Rev. 27(3), 304–311 (2020). [CrossRef]  

27. G. Enemali, R. Zhang, H. Mccann, and C. Liu, “Cost-Effective Quasi-Parallel Sensing Instrumentation for Industrial Chemical Species Tomography,” IEEE Transactions on Industrial Electronics, (DOI: 10.1109/TIE.2021.3063963).

28. Y. Bao, R. Zhang, G. Enemali, Z. Cao, B. Zhou, H. McCann, and C. Liu, “Relative Entropy Regularized TDLAS Tomography for Robust Temperature Imaging,” IEEE Trans. Instrum. Meas. 70, 1 (2020). [CrossRef]  

29. Y. Feng, J. Chang, X. Chen, Q. Zhang, Z. Wang, J. Sun, and Z. Zhang, “Application of TDM and FDM methods in TDLAS based multi-gas detection,” Opt. Quantum Electron. 53(4), 195 (2021). [CrossRef]  

30. C. S. Goldenstein, C. L. Strand, I. A. Schultz, K. Sun, J. B. Jeffries, and R. K. Hanson, “Fitting of calibration-free scanned-wavelength-modulation spectroscopy spectra for determination of gas properties and absorption lineshapes,” Appl. Opt. 53(3), 356–367 (2014). [CrossRef]  

Data availability

Data underlying the results presented in this paper are not publicly available at this time but may be obtained from the authors upon reasonable request.

Cited By

Optica participates in Crossref's Cited-By Linking service. Citing articles from Optica Publishing Group journals and other participating publishers are listed here.

Alert me when this article is cited.


Figures (14)

Fig. 1.
Fig. 1. Dependence of (a) pn(m) and (b) qn(m) on the modulation index m.
Fig. 2.
Fig. 2. Dependence of kp/kq on d when -0.85 ≤ d ≤ 0.85.
Fig. 3.
Fig. 3. Flow chart for calculating gas properties using the proposed algorithm.
Fig. 4.
Fig. 4. Dependence of the condition number of MTM on modulation index m.
Fig. 5.
Fig. 5. The normalized line-shapes at different temperatures for the transitions at around (a) 7185.60 cm-1 and (b) 7444.36 cm-1, respectively.
Fig. 6.
Fig. 6. Simulated WMS raw signal and the extracted even harmonics acquired by the FFT spectral analysis at 1000 K. (a) and (b) show the modulated absorbance signal for the transitions at around 7185.60 cm-1 and 7444.36 cm-1 with the simulated H2O spectral absorption, respectively. (c) and (d) shows the second- to tenth-order even harmonics extracted from (a) and (b), respectively.
Fig. 7.
Fig. 7. Dependence of the objective function F on the modulation index m at 1000 K for the transitions at (a) 7185.60 cm-1 and (b) 7444.36 cm-1.
Fig. 8.
Fig. 8. Relative errors of the retrieved temperature with different harmonic combinations.
Fig. 9.
Fig. 9. Performance evaluation of the proposed algorithm with noise-contaminated data. (a) and (b) show the average relative error and relative standard deviation (RSD) of the temperature retrieved from 100 repetitive simulations.
Fig. 10.
Fig. 10. Schematic diagram of the experimental system.
Fig. 11.
Fig. 11. The experimentally sampled (a) background intensity and (b) transmitted intensity signals at 973 K.
Fig. 12.
Fig. 12. The even-order (2nd to 10th) harmonics demodulated by FFT spectral analysis for the transitions at (a) 7185.60 cm-1 and (b) 7444.36 cm-1.
Fig. 13.
Fig. 13. Comparison of the temperature measured using the proposed algorithm and the thermocouple. (a) and (b) show the absolute measured temperature and the relative error, respectively.
Fig. 14.
Fig. 14. Measured and fitted WMS-2f/1f waveform for the transitions at (a) 7185.60 cm-1 and (b) 7444.36 cm-1 with the residuals shown in above subplots.

Tables (2)

Tables Icon

Table 1. Parameters of the selected transitions at around 7185.60 cm-1 and 7444.36 cm-1. For each transition, the table shows its wavenumbers, line strengths (S(T0)) at T0 = 296 K, pressure broadening-induced half-widths for itself (ξself) and air (ξair), energies of the lower-state (E′′), and coefficients of temperature dependent broadening (nair) [25].

Tables Icon

Table 2. Summary of computational time with different algorithms.

Equations (31)

Equations on this page are rendered with MathJax. Learn more.

v = v 0 + a cos ( ω t ) ,
I t = I ¯ 0 ( 1 + k = 1 i k cos ( k ω t + ψ k ) ) ,
α ( v ) = ln ( I t / I 0 ) = A φ ( v ) ,
φ ( v ) = c L φ L ( v ) + c G φ G ( v ) ,
φ L ( v ) = 2 / ( π λ ) [ 1 + ( ( v v 0 ) / ( λ / 2 ) ) 2 ] 1 ,
φ G ( v ) = 2 / λ ln 2 / π exp [ ln 2 ( ( v v 0 ) / ( λ / 2 ) ) 2 ] ,
c L = 0.6818817 + 0.6129331 d 0.1838439 d 2 0.1156844 d 3 ,
c G = 0.3246017 0.6182531 d + 0.1768139 d 2 + 0.1210944 d 3 ,
λ = 0.5346 λ L + 0.2166 λ L 2 + λ G 2 ,
d = ( λ L λ G ) / ( λ L + λ G ) .
φ ( v ) = 2 π λ [ c L L ( Δ , m ) + c G π ln 2 G ( Δ , m ) ] ,
L ( Δ , m ) = [ 1 + ( Δ + m cos ( ω t ) ) 2 ] 1 ,
G ( Δ , m ) = exp [ ln 2 ( Δ + m cos ( ω t ) ) 2 ] .
L n ( Δ ) = ε n i n 2 m n [ ( 1 i Δ ) 2 + m 2 ( 1 i Δ ) ] n ( 1 i Δ ) 2 + m 2 + c . c . ,
h n L = { ε n m n ( 1 + m 2 1 ) n 1 + m 2 for even  n 0 for odd  n .
h n G = { 2 exp ( m 2 ln 2 / 2 ) I n / 2 ( m 2 ln 2 / 2 ) for even  n 0 for odd  n ,
h n = { 4 A / π λ c L p n ( m ) + 4 A / π λ c G q n ( m ) for even  n 0 for odd  n ,
p n ( m ) = ε n 2 m n ( 1 + m 2 1 ) n 1 + m 2 ,
q n ( m ) = π ln 2 exp ( m 2 ln 2 / 2 ) I n / 2 ( m 2 ln 2 / 2 ) ,
H = k p P ( m ) + k q Q ( m ) ,
H H , = ( E M ( M T M ) 1 M T ) H = 0 ,
m = arg min m | | ( E M ( M T M ) 1 M T ) H | | 2 .
[ k p k q ] = ( M T M ) 1 M T H ,
k p k q = c L c G = 0.6818817 + 0.6129331 d 0.1838439 d 2 0.1156844 d 3 0.3246017 0.6182531 d + 0.1768139 d 2 + 0.1210944 d 3 .
λ = 2 a / m ,
A = π λ k p / ( 4 c L )  or  A = π λ k q / ( 4 c G ) .
H = { k p P ( m ) for 0 .85 <  d     1 , k q Q ( m ) for - 1 d  < - 0 .85,
m = { arg min  m H P ( m ) | | P ( m ) | | 2 for 0 .85 <  d     1 , arg min  m H Q ( m ) | | Q ( m ) | | 2 for - 1 d  < - 0 .85 .
k p = | | H | | 2 | | P ( m ) | | 2 ,
k q = | | H | | 2 | | Q ( m ) | | 2 .
α ( v ) = ln ( I 0 exp ( α ( v ) ) + n o i s e t I 0 + n o i s e 0 ) ,
Select as filters


Select Topics Cancel
© Copyright 2024 | Optica Publishing Group. All rights reserved, including rights for text and data mining and training of artificial technologies or similar technologies.