1. Introduction

Precise and fast measurements of the ¹³C/¹²C ratio in carbon dioxide (CO₂) is important in a number of different areas, such as the medical, biological, environmental and geological sciences. These applications stimulate the development of novel experimental techniques which are capable of high precision and accuracy measurements and short data acquisition time. As far as the measurement accuracy and precision are concerned the instruments measurement precise is the first factor to be considered. The best results can be obtained with isotope ratio mass spectrometry (IRMS), which generally achieves a very high precision in the 0.01‰–0.05‰ range. However, the IRMS cannot differentiate isotopes with almost the same mass, it often requires a complex sample preparation and the precision IRMS is not readily field deployable for in situ measurements. These conditions are very difficult to achieve outside of a well-conditioned laboratory. Development of optical isotope techniques over the last several years has provided scientists a set of techniques for the real-time measurements of the ¹³C/¹²C ratio based on optical absorption which rely on the absorption of infrared laser radiation in vibrational-rotational transitions. These techniques include non-dispersive infrared spectroscopy [1], Fourier transform infrared spectroscopy [2], opto-galvanic spectroscopy [3], photoacoustic absorption spectroscopy [4], cavity ringdown spectroscopy [5], integrated cavity output spectroscopy (ICOS), cavity enhanced absorption spectroscopy (CEAS) [6], direct absorption spectroscopy and wavelength modulation spectroscopy [7]. In these techniques, different isotopic molecular species can be distinguished easily irrespective of their masses by selecting appropriate absorption lines and, the high signal-to-noise ratio (SNR) can be achieved with a high-finesse optical cavity or a multi-pass cell. However, the accuracy and preciseness of the isotope measurement has been affected by the intrinsic noises and interferences from electrical components, background changes and wires as well as environmental changes. The intrinsic noise such as flicker noise, shot noise and detector noise are difficult to remove by hardware system. Many methods, such as least-squares method, background subtraction, the wavelength modulation technique and adaptive filters, have been proposed for this aim to increase measurement system performance. Though these techniques permit one to obtain a high measurement precision that can be similar to that obtained by the IRMS method [8,9], the drawback is a slower response of the system to fast changes in the isotopic ratios. As a result, there is a need for smaller, portable, less labor-intensive, the online monitoring and relatively cheaper instruments with more efficient filtering technology that can operate in the field and make frequent measurements of both the absolute concentration of CO₂ and its isotopic composition. For the specific applications, it is highly desirable to be able to perform real-time measurements with high sensitivity and precision, and to maintain a fast system response. In recent years, digital filtering techniques of wavelet transform (WT) and Kalman filter(KF), which are efficient tools to reduce noises from severely polluted signal without additional hardware, have aroused the interest of researchers and been successfully applied to various spectroscopic techniques for signal processing [10,11], such as Fourier transform infrared spectroscopy [12], laser-induced breakdown spectroscopy [13], Tunable diode laser absorption spectroscopy (TDLAS) [14–17] and wavelength modulation spectroscopy [18],as well as Raman spectroscopy [19].

As superior signal processing techniques, the WT and the KF have been introduced for the isotope measurement [11,15]. In gas isotope measurement, only one of the two techniques is used to suppress the noise. It is still a challenge to get qualitative or quantitative analysis between the WT and the KF, because various noise sources make the denoising process extremely difficult. The comparison of denoising results between different wavelets has not been reported and the comparison between WT and KF have not also been reported.

In this work, a brief description of the TDLAS apparatus for the measurement of δ¹³CO₂ and experiment parameters are introduced in section 2. The signal denoising method and analysis that combines the WT and the Kalman filter are analyzed in section 3. The denoising results of the WT and the KF are quantitatively compared and the KF is selected an efficient suppressing noise tool for the measurement of δ¹³CO₂ as presented in section 4.

2. Experimental setup and noise analysis

The TDLAS is used for measurement of δ¹³CO₂ and the scheme of the experimental setup is depicted in Fig. 1. The system consists of three main parts, named as the optical system, detector system and electronic system. The laser source was a cw room temperature near-infrared distributed feedback (DFB) diode laser operating at 2.008 µm and, it can be tuned in the range from 4976 to 4984 cm⁻¹ based on our analysis of absorption lines intensities at a fixed temperature of 24 °C. In this spectral region, two relatively strong absorption lines of the ¹²CO₂ vibrational bands ν₁ + 3ν₂ + ν₃−ν₂ (R17 13¹1 ← 01¹0, ν = 4978.204746 cm⁻¹) and ¹³CO₂ vibrational band ν₁ + 2ν₂ + ν₃ (P16 12°1 ← 00°0, ν = 4978.022037 cm⁻¹) are selected for avoiding H₂O disturbance. The laser beam is directed to a compact dense-pattern multiple-pass cell with an optical path length of 29 m. The laser spatial propagating characteristics and the spot distribution in the multi-pass cell are shown in embedded graph of Fig. 1. From Fig. 1 it can be found that the optical interference fringes (etalons) can be neglected. Then the signal is detected by PDA10DT(-EC) (Thorlabs). A sample temperature stabilization to better than ± 0.01 K is, thus, required ensuring a measurement precision <1‰. The temperature of the multiple-pass cell is actively controlled and can be maintained constant within ± 0.01°C at the set temperature. The cell temperature is monitored by three calibrated thermistors (PT100) (with an accuracy of 0.03 °C and a precision of 0.01 °C). The 1000 ppm calibrated CO₂ is measured in the following experiments.

 

Fig. 1 Scheme of the experimental setup.

Download Full Size | PDF

The intrinsic sources of noise in TDLAS system can be mainly divided into photon shot noise, the dark current and the flicker noise, as well as Johnson noise (thermal noise). Photon shot noise, which is an inherent noise source and results from the random nature of photon emission, is present in the optical system. Although the integration time and quantum efficiency of the detector influence photon shot noise, it is impossible to avoid for its intrinsically random nature. The dark current is generated by the random generation of electrons and holes in the depletion region of the photodetector in the absence of light. It can be minimized by cooling the photodetector. The flicker noise can be negligible for the attenuation of the flicker noise in the most of avalanche diodes. Johnson noise (thermal noise) is generated by the random thermal motion of electrons and it happens regardless of any applied voltage. Johnson noise is approximately white noise and its power spectral density is nearly equal throughout the frequency spectrum and can't also be avoided. These noise can hardly been removed by increasing the length of the absorption cell. In practical applications, it is necessary to study simple and effective ways for these noise suppression. The digital filter technique of the WT and the KF can be used to suppress the noise without increasing the volume of the system. The δ¹³CO₂ was measured by the linear regression method [20], and the application of the single absorption cell can eliminate the influence of optical path difference on the measurement results. When the δ¹³CO₂ are measured, experiment data collection, background subtraction, window selection, linear regression, and real-time signal processing with digital filtering techniques of the WT and the KF are integrated in the same program of LabVIEW software, which implements the real time, continuous, and coordinated operation.

3. Analysis and results

3.1 Wavelet transform denoising

Wavelets were introduced in early 90s with numerous applications in electroencephalogram [21]. In recent years, the WT, which shows its capacity in denoising and the enhancement of SNR, has been successfully applied to various spectroscopic techniques for signal processing [12–19]. The theory of the WT denoising has been well established in the spectrum measurement [22,23]. The wavelet threshold denoising method is proposed in the trace gas isotope measurement and can be divided into three steps as following, wavelet basis select, wavelet function decomposition, optimal threshold and the final estimated denoising results.

In the process of wavelet denoising, it itself does not produce a reduced version of the original signal. It is used by coupling with a feature selection strategy and a subset of wavelet coefficients with valuable information are picked out. We always hope that the selected wavelet basis function can meet the properties of 1) orthogonality, 2) high vanishing moment and 3) the compactness. However, the selected wavelet basis function can hardly satisfy above the three properties. In the practice application, the appropriate wavelet basis function will be selected according to the specific measurement requirements. Sometimes, it is necessary to construct a orthogonal wavelet basis, however, the orthogonal wavelet bases can’t simultaneously improve the resolution in time domain and frequency domain. With the increase of decompositon level, the frequency window becomes wider, which results in a low resolution in frequency domain. In order to improve the resolution in frequency domain, the widen frequency window needs to segment. The wavelet packets are introduced to satisfy this purpose. The signal can be decomposed more precisely by use of the wavelet packets. In the measurement of δ¹³CO_2, all wavelet-based approaches and KF are implemented using LabView 8.2 under Microsoft Windows environment, which includes wavelet packets. The choice of wavelet basis function directly affects the wavelet packet transformation processing result to signal. A variety of wavelet basis function can be used for wavelet packet analysis of the δ¹³CO₂ measurement signal and the optimal wavelet basis function is selected. Because the WT has the characteristic of multi-resolution analysis, the signal can be decomposed into a series of sub bands with different frequency, and the high frequency signal is separated from the low frequency signal. Different wavelet basis functions are used to analyze the same absorption signal and the denoising results will be obtained. For exhaustive comparative analysis, we empirically investigate the performances of the most five common wavelet basis function namely, Harr, Daubechies, Symlets wavelet, Coiflets wavelet, and Biorthogonal Wavele and the results are shown in Fig. 2. Through the comparison of the various wavelet basis functions and achieving the best results according to the smoothness of the isotope spectrum measurement, the optimal wavelet basis function of Haar will be used to suppress the noises.

 

Fig. 2 Comparison of the Original Signal and the resulted signal of the different wavelet basis function.

Download Full Size | PDF

The wavelet decomposition level (DL) is another important process after the optimal wavelet basis function is determined. In order to select the optimal DL and effectively better noise reduction of the δ¹³CO₂ measurement, one objective is to put forward a scheme of choosing optimal DL. In theory, The maximum decomposition scale of wavelet basis function is $J=\left|{\log }_{2}N\right|$, represents the integer and the value of J is 3—7. Wavelet decomposition of Harr wavelet basis function is performed and the results are shown in Fig. 3. In theory, when J is larger, the difference between the signal and noise is more obvious, the noise is easy to suppress. From the Fig. 3(a) it can be seen that with the increase of wavelet DL, the high frequency noise of the spectrum is removed step by step and, the smoothness of the spectrum is getting better and better. It also can be seen from the Fig. 3(a) when the DL is too small, the noise suppressing results of measured signal is not obvious so the result of noise reduction is not ideal. When DL is too excessive, the signal distortion and the loss of useful information will be caused and the measurement precise will decrease as shown in Fig. 3(b). From Fig. 3(b), it can be found that the effect on the weaker (¹³CO₂) is stronger than the stronger (¹²CO₂) when the DL exceeds the 5th level. Which results in decrease of the measurement precise of the δ¹³CO₂. Therefore the noises are removed basis on the optimal DL of the fifth level.

 

Fig. 3 Comparison of different DL denoising results based on Haar wavelet basis: (a) the spectrum smoothness after denoising; (b) the spectrum change of ¹²CO₂ and ¹³CO₂ before and after denoising, the left spectrum is ¹²CO₂ absorption spectrum in the range from 471—517 and the inset of the upper panel is enlarged signals in the range from 485—500, the right spectrum is ¹³CO₂ absorption spectrum in the range from 1260—1311.

Download Full Size | PDF

After the wavelet basis function and the optimal DL are determined, the next step in wavelet based denoise is wavelet threshold and the threhold function selection. The concept of denoising with wavelet threshold, based on the WT was brought out by Donoho and Johnstone in 1994. In the WT denoising methods, the wavelet threshold denoising (WTD) method is the simplest in realization and it requires the least computational time. Thus the WTD has been used widely. In gas isotope measurement the threshold selection has directly relations with the denoising results. The different denoising results can be obtained with the different wavelet threshold when the same wavelet function is used to suppress the noise. Partial wavelet coefficients can't be set zero when the threshold is undersize so parts of noises are retained and some useful signals will be taked off when the threshold is iusto major so parts of useful signals are lost. Therefore, how to select effectively threshold to avoid lose useful signal is problem worth to research. In the same denoising condition, typical threshold of VisuShrink (VISU), Sureshrink (SURE), Hybrid and Minimax are used in the experiment and, the denoise results are shown in Fig. 4, and singular point is indicated in the dashed frame. From the Fig. 4 it can be found that the denoise effect are less effective for the spectral singular point by use of the SURE threshold and the Minimax threshold. By comparing the spectral smoothness we can see that the denoising effect by use of the VISU threshold is better than that of Hybrid threshold. Therefore, in the δ¹³CO₂ measurement the VISU threshold is used to suppress the noise.

 

Fig. 4 The comparison of the denoising results with different thresholds.

Download Full Size | PDF

Wavelet threshold faces the question of the threshold function selection in its application. The value larger than threshold are regarded as signal, and the smaller ones are regarded as noise. In the measurement of gas isotopes, the global threshold function and the local threshold function can be used to estimate the denoising results. The global threshold function denoising has the same the wavelet coefficient for all the DL or the same DL. The local threshold function is more flexible than the global threshold in noise reduction, which can select the appropriate threshold based on the current wavelet coefficient. In the case of the spectrum denoising, the local threshold function is selected and the measurement results are shown in Fig. 5. From Fig. 5 it can be found that the denoising effect of the local threshold function is better than that of the global threshold function.

 

Fig. 5 Estimation with the global threshold and the local threshold. The inset of the upper panel is the enlarged signals in the range from 865 to 872.

Download Full Size | PDF

3.2 Kalman filter

The intrinsic noise of the spectrometer can be reduced by selecting stronger absorption lines and using a multi-pass cell or a high-finesse optical cavity. Another possibility is to selected the optimal averaging time which can be determined by an Allan variance analysis [24]. The Allan variance defines the optimal averaging time for achieving high precision. Once this averaging time has been optimized, the high measurement precision is achieved. Further measurement precision improvements can be achieved through the use of the KF technique.

The KF was first developed in 1960 and applied to aerospace navigation problems [25] and, applied to tunable diode-laser-based gas sensors was studied by Riris et al. in 1994 [26]. In recent years the KF has been an optimal state estimation process applied to a dynamic system and has been applied to real-time trace gas concentration measurements [27] and the measurement of H₂O isotopologue ratios [11]. Compared with the simple multi-point averaging method, the advantages of the KF technique can be described as:

Then the KF technique is applied to fast and high-precision isotope ratio determinations. The the detailed theoretical analysis of the KF technique that is applied to absorption spectrum measurement can be found in [11,26,27]. Basis on the linear stochastic difference model, the KF can be two steps of prediction and correction. The brief description of the KF model for isotope ratios measurements is given as following.

The true isotope ratio value δ_k+1 at time k + 1 is evolved from the value δ_k given at time k according to

Before experiment a priori values of ${\stackrel{⌢}{\delta }}_k^{-}$ for δ and the error variance $e_k^{-}$ are selected. The selection of these quantities are not critical, since the filter will converge to an appropriate value. However, they should be chosen to be within the dynamic range of the expected isotope ratio and errors to speed convergence. The next step involves taking a measurement update from the sensor by first computing the Kalman gain. Subsequently, the sensor performs a isotope ratio measurement to obtain δ. This is followed by calculating an a posteriori δ estimate, and the final step is to obtain an a posteriori error variance estimate. Once a posteriori estimates of the δ and variance are obtained, a time update is processed, if a control input is modeled. This step projects the isotope ratio and variance estimates from time step k to step k + 1. This recursive nature is one of the appealing features of the KF and, the KF is realized by a laptop PC running LabView 8.2 software.

We firstly use conventional averaging of n measured δ values, and the optimum number n will be determined from an Allan variance analysis. The δ¹³CO₂ was measured within 5000 s and the Allan variance plotted as a function of the averaging time leads to the Allan plot as shown in Fig. 6. From the Fig. 6 we can see that the optimal averaging time of the system is about 90 s. Then the KF is applied to measure the δ¹³CO₂ with the sampling period of 1 s shown in Fig. 7(b) and the measurement precise of the δ¹³CO₂ at the optimal averaging time is shown in Fig. 7(a). From the Fig. 7 it can be seen that the precise of the δ¹³CO₂ at the optimal averaging time is higher than that of the KF, however, this improvement, obtained by post processing of the data, comes at the price of a slow system response and the detail change of the δ¹³CO₂.

 

Fig. 6 Upper panel, raw measurements of δ¹³CO₂ with 1 s averaging time, the Allan variances in the lower panel show an optimal averaging time of about 90 s for the present laser system.

Download Full Size | PDF

 

Fig. 7 The δ¹³CO₂ measurements showing a comparison between the optimal averaging solutions and the KF solutions.

Download Full Size | PDF

In order to quantitatively compare the denoising results of the WT and the KF for the δ¹³CO₂, the calibrated 1000 ppm CO₂ is measured in the same experimental condition, and the results are shown in Fig. 8. The top panel is the original signal of simple averaging, the middle panel is the WT denoising signal and the bottom panel is the KF denoising signal. From the Fig. 8 it can be found that the original achieved the worst results and, the WT and the KF get similar results. The WT and the KF offer performance that is superior to simple averaging. While the WT results in superior denoising, this comes at the price of the measurement time caused by the selection of the wavelet basis function, the optimal DL and the wavelet threshold.

 

Fig. 8 The comparison of δ¹³CO₂ measurement between the original results, WT results and KF results. The top panel is the original signal, the middle panel is the WT denoising signal and the bottom panel is the KF denoising signal.

Download Full Size | PDF

4. Conclusion

In this paper, the WT and the KF are successfully used in suppressing noise for δ¹³CO₂, and the experiments shows that they have been effective data analysis tool for signal processing and can get satisfactory results. According to the denoising results between the WT and the KF, they give 1σ similar precision of 0.0666‰ and 0.0619‰ for δ¹³CO₂ shown in Fig. 8, respectively.

Though the WT denoising method is a powerful denoising tool for trace gas isotope ratio measurement, its effectiveness is influenced by the selection of the optimal wavelet basis function, the optimal DL and the wavelet threshold. These factors have to be taken into account in suppressing noise with the WT. It takes long time to use wavelet to suppress the noise.

The KF can be applied in two steps of prediction and correction. The denoising process of the KF is relatively simple. The KF is an optimal recursive data processing algorithm and used to estimate the performance of a system from measurements which contain random noise. In the authors' opinion, the KF can be chosen by using the proposed method for real-time, online measurements of trace gas isotope ratios.