Login| Sign Up| Help| Contact|

Patent Searching and Data


Title:
MULTICOMPONENT ANALYSIS OF FT-IR SPECTRA
Document Type and Number:
WIPO Patent Application WO/1993/025890
Kind Code:
A1
Abstract:
A procedure for analyzing a multicomponent FT-IR spectrum of an unknown mixture of gases, wherein coefficients xj of Eq. (1) is determined and the partial pressures of the pure gases in the mixture are obtained by multiplying the measuring pressures of the corresponding library spectra by their coefficients in the linear combination; in which procedure a solution of the optimization proble of Eq. (4) is obtained through fulfilling the coefficient vector (7); the inner product <a|b> is defined according to Eq. (8) and the Eq. (7) is settled according to Eq. (9).

Inventors:
KAUPPINEN JYRKI (FI)
SAARINEN PEKKA (FI)
Application Number:
PCT/FI1993/000252
Publication Date:
December 23, 1993
Filing Date:
June 09, 1993
Export Citation:
Click for automatic bibliography generation   Help
Assignee:
TEMET INSTR OY (FI)
KAUPPINEN JYRKI (FI)
SAARINEN PEKKA (FI)
International Classes:
G01J3/45; G01N21/3504; (IPC1-7): G01N21/35
Domestic Patent References:
WO1991019183A11991-12-12
Foreign References:
US4928015A1990-05-22
EP0305090A21989-03-01
US4500207A1985-02-19
US5113073A1992-05-12
Download PDF:
Claims:
Claims
1. A procedure for analysing a multicomponent FTIR spectrum of an unknown mixture of gases by an analyzer in which a set of FTIR library spectra of pure molecular gases, measured in known pressures, are stored as vectors k having N components and by which analyzer inner products for all possible combinations of two vectors k are calcu¬ lated and stored, the procedure comprising the steps of choosing of M library vectors k for the analysis of the unknown mixture of gases forming thus a basis matrix K, forming an orthogonal basis matrix K' by deforming the vectors k chosen for the matrix K to orthogonal vec¬ tors k' so that inner products of these vectors k' with each other but not with themselves are zero and calculat¬ ing and storing the inner products of these orthogonal vectors k' with themselves, forming and storing a basis transformation matrix Q so that K = K'Q, forming and storing an inverse matrix Q"1 of the matrix Q, measuring an FTIR spectrum s for the unknown mix ture of gases as a vector having N components, calculating an optimal coefficient vector x' for the spectrum s in the orthogonal basis by using for its components x'j the equation XJ (k 'J \ k 'J) calculating a coefficient vector x so that x = Q"1 x', and calculating partial pressures of the pure gases in the mixture of gases to be analysed by multiplying the measuring pressures of the corresponding library spectra by their coefficients xi thus forming an analysis result for the spectrum s.
2. A procedure according to claim 1 further comp¬ rising comparing the measured spectrum s to its analysis result and in the case of mismatch adding a new library vector k M+l to the basis matrix K by calculating a new orthogonal vector k' M+l, calculating and storing an inner product of the new orthogonal vector k'M+1 with itself, adding a new column to the transformation matrix Q and storing this new matrix, forming and storing a new inverse matrix Q"1 on the basis of the new matrix Q by adding a ne column to the inverse matrix Q"1, calculating a new component x'M+1 for the orthogonal coefficient vector x' of the spectrum s to be analyzed b the equation χ _ (k 'M+1 \ s) M+1 ik' l k""1) and adding it to the vector x' , and calculating a new coefficient vector x so that x Q"1 x'.
3. A procedure according to claim 1 or 2 furthe comprising evaluating a deviation σs for the noise of th spectrum s, calculating variances for the coefficients x^ of th vector x by using the equation and calculating error limits vi (p) for the coefficient x. by the equation Vj ~ JVa ljc i (i J2) . in which p is the probability by which the error for the library spectrum j is between the limits (Vj, Vj) and Φ1 is the inverse of the distribution function of the one dimensional standardized normal distribution.
4. A procedure according to claims 2 and 3 further comprising calculating a variance for the new coefficient XH+1 by the equation Vaτ (x^. ) = *+1 N ( k'M^ \ k'M+1) and calculating values to be added to the variances of the coefficients x.., j = 1,...,M by the equation N (it'**1!*'"*1) ' adding said values to the variances of the coeffi¬ cients xjf j = 1,...,M and calculating new error limits v., (p) for the coeffi¬ cients Xj by the equation in which p is the probability by which the error for the library spectrum j is between the limits (vjf v^ ) and Φ"1 is the inverse of the distribution function of the one dimensional standardized normal distribution.
Description:
MULTICOMPONENT ANALYSIS OF FT-IR SPECTRA

Background Of The Invention

The present invention concerns the procedure for analyzing a FT-IR-spectrum of an unknown mixture of gases,about which we do not know the constituent gases, not to mention their partial pressures. Instead we know a large set of library spectra of pure molecular gases, measured in known pressures and with the same interfe¬ rometer as the unknown sample.

1. Field Of The Invention

In FT-IR-spectroscopy in multicomponent ana¬ lysis we have the IR-spectrum of an unknown mixture of gases, about which we do not know the constituent ga¬ ses, not to mention their partial pressures. Instead we know a large set of library spectra of pure molecular gases, measured in known pressures and with the same interferometer as the unknown sample. By using these pure spect ra we ought to calculate the partial pressu¬ res of the pure gases in the mixture, with error li¬ mits. The errors in the obtained values arise from the measurement noise in the spectra. Because the calcula¬ tion of such partial pressures that best explain the mixture spectrum is a relatively simple task, we shall consider it only briefly and concentrate on the calcu¬ lation of their error limits. We shall also consider the optimal choice of the resolution providing as small error limits as possible and the application of the non-negativity constraint for the partial pressures. In addition to gas spectra, all the calculations apply to spectra of non-interacting liquids as well.

2. Description Of The Prior Art

Let us denote by s the measured mixture spect¬ rum to be analyzed and let {K j |j = 1,...,M} be the set of library spectra with which we are trying to explain s. We assume that all the spectra are linearized, i.e. that on the y-axis there is always the negative loga¬ rithm of the transmittance (except possibly some cons¬ tant coefficient) . This means that we are using the absorbance scale. Then, according to Beer's law, s is some linear combination of the spectra K j , or

Our task now is to determine such coefficients by which s is explained as well as possible. After that the partial pressures of the pure gases in the mixture are obtained by multiplying the measuring pressures of the corresponding library spectra by their coefficients in the linear combination.

Because all the spectra are only known at equally spaced tabelling points υ,,...,υ N , s and K ~ can be treated as vectors, and we can state Eq. 1 in vector form as

This equation can be further simplified by collecting aallll KK jj --vveeccttoorrss ttooggeetther into a matrix K, which gives us the matrix equation

Kx = s. (3)

where

and x is the coefficient (column) vector containing the coefficients κ. of the pure spectra.

Summary Of The Invention

Now all the spectra that can be exactly ex¬ plained by the used set of M pure spectra are expressi¬ ble in the form Kx. However, even if the set {K^lj = 1,...,M} included all the components in the mixture, there always exists some noise in the spectra so that s never lies exactly in the vector space generated by the vectors K 3 . Thus we have to content ourselves with such a linear combination Kx of the pure spectra that is as near as possible to the vector s. This means that we have to minimize the norm of the residual vector s-Kx, which remains unexplained. Thus we can write our prob¬ lem in the form of an optimization problem as

II s - Kx II = min! . (4)

The optimization parameters are then the components of the coefficient vector x.

The optimal value of the coefficient vector depends on which vector norm we are using. Let us now define this norm from the condition

This is the usual square sum vector norm except for the

constant coefficient 1/N. The purpose of this coeffi¬ cient is to make vector norms as independent as possi¬ ble of the number of the data N (and thus on the reso¬ lution) . It has no effect on the solution of the opti¬ mization problem 4. The solution of problem 4 is now obtained by setting the gradient of the residual vector equal to zero or

V || s - Kx || = 0. (6)

a very simple calculation shows that this condition is fulfilled by the coefficient vector

x opt = (K τ K) _1 K τ s (7)

(provided, of course, that N > M. ) This is known in literature as the pseudoinverse solution of the matrix equation Kx = s (see e.g. [1]). Because of the choice -see Eq. 5- of the norm, this can also be called the least squares solution. The matrix K T K is invertible, if and only if the vectors K 3 are linearly independent. If they are not, the best solution is not unique. If the error limits are not needed, we suggest the procedure described below for the application of Eq. 7. We now define the inner product (a|b) of two vectors a and b from

The coefficient 1/N again ensures that the inner pro¬ duct is as independent as possible of the number of the data. Likewise the usual relationship || v ||= (v|v) 12 between norm and inner product holds true. Now we see that equation 7 can be written in the form

where

and

Y j = <K j |s) . (11)

Thus the matrix A contains the inner products of the used library spectra with each other and may be called the inner product matrix. If we are not interested in the error limits of the coefficients, we can calculate beforehand all possible inner products into a library. Then, in a particular analysis, we need only pick out those inner products corresponding to the set of libra¬ ry spectra to be used in this analysis to build up the matrix A, calculate the vector y and solve the matrix equation Ax = y. If this is done by using the Gauss elimination method (see e.g. [2, p. 179]), there is no need to invert A, and the analysis is a very rapid procedure.

All noise in the spectra generate errors in the optimal coefficient vector of Eq. 9. In this work we shall only deal with the noise in the mixture spect- rum s. The noise in the library spectra K j could appro¬ ximately be taken into account, too, but we shall omit it here for two reasons. For the first the treatise would be too long and only approximate. For the second, in practical applications it is sensible to register the library spectra with much higher accuracy (i.e. by co-adding a much larger amount of individual runs) than the mixture spectrum, which is obtained with a rapid measurement. Thus, practically seen, all the noise is in the mixture spectra and the library spectra can be considered as noiseless.

Above we have been using the library spectra 3 as the basis vectors, with which all measurements were

explained in the form Σ M 3-1κ3.K 3 or Kx. However, when calcu- lating the error limits, it is advantageous to use an orthogonal basis, which generates the same vector space as the original library spectra. We construct these orthogonal basis vectors by the Gram-Schmidt ortho- gonalization method (see some textbook of linear algeb¬ ra, e.g. [2, p. 138]). This means that the new basis vectors {K' j } are defined by the recursive group of equations

Thus the old basis vectors are obtained from the new ones by

K^K' 1

Next we define a transformation matrix Q as

(Thus Q is an upper triangular matrix.) Now the trans¬ formation of the basis can be expressed very compactly as

K = K'Q. (15)

The transformation of coefficient vectors between these two bases is governed by the equation

x = Q^x' . (16)

Thus we shall need the inverse of Q. After some tedious calculations we discover that it is an upper triangular matrix given by

0, . i > j

The elements with j>i have to be calculated in the order i = M - 1,...,1.

The determination of the best coefficient vector is very simple in the basis K' . Because this basis is orthogonal, (K !i |K ,:1 ) = 5.. and the inner pro¬ duct matrix A' (see (Eq. 10)) is diagonal.

A' = diag( || K' x" 1 w K' ' (18)

Thus the best coefficient vector x 1 opt is given by (see

or

χl _ <K"|s> ( 20 )

'opt, ,, κ [J2

Vector x' ot gives the optimal coefficients of the orthogonal basis vectors K ,j . The corresponding coefficient vector x of the real pure spectra K j is then Q ~ x' . As is intuitively very clear, this equals the optimal coefficient vector x in the original basis. This fact can also be proved formally as follows:

x = [(ζT) x ICKQ- (Q ) Ks = = Q(K T K) _1 K T S = QX opt .

This means that

Next we consider the effect of measurement noise on the optimal coefficient vector of Eq. 9 or 21. For this purpose we divide s into two separate parts as

s = s° + s e , (22)

where s° is the correct noiseless mixture spectrum and s e includes the noise. Because the dependence between s and xop_t is linear,' the coefficient vector xop _t can also be divided into the correct coefficient vector xo°p„t and the error vector x e , which obey

xV opt = (K^Γ ' VS 0 = Q (K ' K ' ) K ' T s 0 ( 23a ) and χ e opt = ( v K T K ) ' _1 K T s e = ( 24a )

or equivalently

and

Thus the errors of the coefficients depend linearly on the measurements. The components of the noise vector s e

are normally distributed with zero mean. Even if the noise data of an individual measurement were not nor¬ mally distributed for some reason or other, according to the central limit theorem of the probability calcu- lus the sum of several noise data is always normally distributed, and in practice we always co-add several individual runs to obtain one spectrum. Thus, if we designate their standard deviation by σ , we can write

s* ÷ N(0,σ ) . (25)

(This expression tells that s* obeys the normal distri¬ bution with a mean 0 and a variance as 2 . Note that N stands for the normal distribution, while N is the number of the data.) Because we know the distribution of the noise data, we can now calculate the distribu¬ tions of the errors of the coefficients by using Eq. 24. This is what we are going to do next.

First we state without proof a well-known result of the probability calculus concerning the nor¬ mal distribution.

Let z ÷ N( μ i , σ 2 i ) be independent and a ± ε R. Then

α - z i ÷ N (∑ « .-/*. > ∑ α χ ) (26)

This means that a linear combination of independent normally distributed random variables is normally dist¬ ributed, too. By using this result we get for the com- ponents of the vector y' e , which appears in the expres¬ sion for x e Op„t in Eq -. 24b:

y" = cκ")*« * + N(o. i _ N(O. # ιι v- f ) < 27 ,

Now we need another result of the probability calculus, stating 'that two normally distributed random variables z and z 2 are independent if and only if the correlation

between them is zero. The correlation is defined by

where

cov( Zl ,z 2 ) = E[( Zl - E( Zl ))(z 2 - E(z 2 ))]

and E is the expectation operator. Because the expecta¬ tion values Efy'*) are zeroes, two components y 1 ^ and y* are independent random variables if and only if E(y' Now

W j V'D = f 2 E[(.K' i s e )(K' k s < )]

,

for E(s ® s ® ) = δ u E[(s*) 2 ] and all the random variables s* are identically distributed. The independence of y 1 ' and y'* now follows from the fact that the basis K' is orthogonal. It further follows from Eqs. 19, 26 and 27 that

' t,: N(0,

N|| K" ||= ) (28)

and that the components of '* t are independent. Thus we may once more use result 26 to obtain

x X < ot,j -~ -=( WQ- 1 x x 'p Λt . •χ I ( Iθ υ,d N 2γ_, | i | Q κ ~7 , l ι l)2 )i • < 29)

If we use the non-orthogonal basis K, the left-hand equality of Eq. 24 gives

M

< 30 > where P is the pseudoinverse matrix

BST SHEET

P = (K T K) 1 -I

(31)

Now that we know the distributions of the errors of the coefficients, we are also able to give error limits for the coefficients. Because these errors are normally distributed random variables, we cannot give any upper limits for them. Instead we can calcula¬ te such error limits for the correct coefficients to fall with a given probability within them. Let us now assume that we want such an error limit υ , that the coefficient κ o t belongs with a probability p in the interval [ 1 κ o°pt,j, - u j j; 0 opt,j. + υ j J . This means that

P ( κ αpt,j [ -«1 » a ] > -

If we denote Φ the distribution function of the stan¬ dard normal distribution (N(0,1)), we get.

from which

For example, the 50% error limit for the coefficient of library spectrum number j is given by the expression

The variance of κ o e pt, j α is obtained from result 29 Or 30 as

2 M (Q-^ ~ M

Varfx' ,) = Z__ y ΣllLL _ σ 2 - P 2

- °P*, j J N N Z £__< |-| ~ κ - K ,,li |n 2 - σ > ___ Cl

/=ι (33)

It might first appear that the original basis K would be the most sensible choice to be used in the analysis, because the coefficients and their variances are then obtained from the very simple results 9 and 30. However, every time we add a new library spectrum into the analysis, the inner product matrix A is chan¬ ged, and we have to re-calculate its inverse, which is needed at least in Eq. 31, and in inverting a general matrix we do not have a simple formula like Eq. 17. In addition we have to re-calculate the matrix product A _1 K T to obtain P. The previous value of x o t also becomes useless. When using the orthogonal basis K 1 , the coef¬ ficients and their variances are given by Eqs. 20, 21 and 29. Let us now have a closer look at what happens when we have to add a new library spectrum to the set of the library spectra used in the analysis and the basis K' is in use.

First of all we must calculate a new ortho- gonal spectrum K' M+1 from group of equations 12. While we calculate the coefficients appearing in the expression K ,M+1 in group of equations 12, we also obtain the ele¬ ments of the column M + 1, which must be added to the transformation matrix Q. (The new row elements Q u , ,, j > M + 1 are zeroes according to Eq. 14.) This new co¬ lumn means that the inverse transformation matrix ζf is changed, too. However, we see from Eq. 17 that an ar¬ bitrary element of Q "1 depends only on the previously calculated elements of Q '1 in the same column and on the elements of Q in the same and previous columns. Thus the previous elements of Q _1 are not changed and we need only calculate M new elements θ M "1 . f M,+1,,...,Q " 1 1 ,M+ 1, which are

given by Eq . 17 as

M+l

Qi,M+ι = ~ 2^ QUQI,M+I (34) l=i+l

(Because Q .-1 is an upper triangular matrix, all the elements on the new row Q M ' + 1 ,l are either zero or one.)'

If we are still analyzing the same measurement s as before updating Q '1 , the following corrections remain to be made to the previous calculations: 1) We have to add a new comp c onent κ' opt,M„+„l into the vec-

VZ M+l ) 2

N || K' M+1 || 2 '

Because the squares of norms || K' j || 2 are nee¬ ded repeatedly, it is sensible to store them in a vec¬ tor, after they have been calculated in group of equ¬ ations 12. The use of double precision real numbers is highly recommendable in computer programs.

Brief Description Of The Drawings

The invention is described in the following detailed in the aid of examples, referring to the draw- ings attached wherein

Fig 1 is a mixture spectrum s, the linear combination Kx o of the library spectra that best explains s and the residual spectrum S-KX Q , which then remains unexplained (appears magnified in the figure); Fig 2 discloses an illustration of what happens when the used library is deficient; the same mixture spect¬ rum has been analyzed as in Fig. 1, but now the spect-

rum of 2-butanone has been removed from the used set of library spectra;

Fig. 3 discloses a situation if the interferometer is properly aligned, but the radiation source is circular instead of an ideal point source, the signal in every individual wavenumber υ 0 is spread uniformly of the interval from υ„ 0( v l-αm 2 ax/2)' to υ„0;'

Fig. 4 discloses a situation when the observed line shape e° is always the convolution of the real line shape e, the boxcar broadening function w Ω described in Fig. 3 and the Fourier transform w L of the boxcar inter- ferogram truncation function;

Fig. 5 discloses a situation when for every component κ of the optimal non-negative coefficient vector x opt either one of the following conditions holds:

(a) 5 i = 0 and x j > 0 or (b) ' S IL > 0 and =0.

Description Of The Referred Embodiment

In this chapter we shall demonstrate the method by analyzing the spectrum of a simple mixture of gases. The coefficients are calculated by using Eqs. 20 and 21 and their error limits by Eqs. 29 and 32. We shall use 50 % error limits here, because they give a very good idea of the orders of magnitude of the er¬ rors. In Fig. 1 we see the mixture spectrum s to be analyzed by using a set of 12 library spectra. In the table 1 we see the components of x with error limits. The three first components are the background spectra. These are all the spectra or functions, which are pres¬ ent also in the background measurement or"which can arise from an error data in the interferogram. They can have negative coefficients as well as positive. The mixture spectrum s and the library spectra K 3 are measu¬ red with identical arrangements, but extra artificial noise has been added to s in order to simulate poorer

accuracy, which is assumed in this work. In addition we see the linear combination KX Q t of these library spectra that best explains s. Also shown is the remainder spec¬ trum s-Kx o t , which then remains unexplained. As we may see, the remainder spectrum consists of pure white noise. This indicates that the analysis has been suc¬ cessful and that the 12 library spectra used suffice to explain the measurement. In the table 1 we can see the result of the analysis or the optimal coefficient vec¬ tor xopt.

The optimal coefficient vector

back- ground

The partial pressure of the library gas number j is now obtained by multiplying its measuring pressure by its coefficient κ o t . Those components, which do not exist in the mixture, have small positive or negative coefficients with the same order of magnitude as their error limits. If we add some background spectra to the analysis, the optimal removal of the background of the .spectrum is also performed automatically. In fact, the three first spectra (water, carbon dioxide and the

constant function, which has the value of one every¬ where) in the table are background spectra. Because we are not interested in their exact coefficients, but merely want to get rid of them, they need not even be pure, but can contain each other. The library spectra may contain some carbon dioxide and water, too, with no effect on the values of their coefficients. The only effect is that the coefficients and error limits of the background spectra are not reliable. The coefficients of the background spectra can have negative values as well as positive. It can actually happen, for instance, that there is more carbon dioxide present during the measurement of the background interferogram than there was during the measurement of the sample interferogram. In Fig. 2 we see what happens if the set of library spectra used is not sufficient for explaining the compound spectrum, in the table 2 we can see the results of the analysis or the optimal coefficient vector κop.t. The optimal coefficient vector

The same mixture spectrum is used as in Fig. 1, but the spectrum of 2-butanon, which is an important constituent of the mixture, is not included in the analysis. Now the minimized residual spectrum is not pure noise and the coefficients of the remaining li¬ brary spectra have changed, too. As we can detect, now

the minimized residual spectrum s-Kx o t is no longer white noise but has a distinctive structure and, which is very important, it has some structure on the same wavenumbers as the missing spectrum has spectral lines. Thus it is possible to deduce from the remainder spect¬ rum what kind of spectra we ought to add to the ana¬ lysis. Because the portion of the missing spectrum has to be explained with the remaining spectra as well as possible, their coefficients are distorted, too, and the error limits are no longer reliable. Thus we must always add new library spectra to the analysis until there is no structure left in the remainder spectrum. If there are different pressure broadenings in s and in the library, it might be helpful to use a couple of library spectra for a single compound. However, a bet¬ ter procedure would be to decrease the resolution so that all lines become sine-function shaped (see the next chapter) .

In Fig. 4 function w L is a sine-function, and the interval between its two successive zeroes equals the data interval Δυ in the spectral domain, which in turn equals l/(2NΔκ). The FWHH of w L is approximately the same. The width W of w Ω is directly proportional to υ Q (the location of the spectral line) and the area of the radiation source. Its height H is inversely propor¬ tional to υ Q and directly proportional to the surface brightness of the source. In the optimal situation the FWHH's of the three right-hand side curves are approxi¬ mately equal. The total area of e° is the product of the areas of e, w Ω , and w L . The area of w L is always 1 and can be regarded as dimensionless.

After we have tixed the wavanumber range we want to deal with, the sampling interval in the inter¬ ferogram is also fixed according to the Nyqvist samp- ling theorem. However, it still remains to choose the length of the registered interferogram or the amplitude of the mirror movement, which in turn is determined by

the used number of the data N. Note that we denote by N the number of the data of the single-sided interfero¬ gram. The corresponding number used by the Fast Fourier Transform algorithm is then 2N. We are now going to examine how the number of the data ought to be selected in order to minimize the error limits. As we can see from Eq. 32, the error limits are directly proportional to the standard deviations of the coefficients, which means the square roots of their variances. From Eq. 29 or 33 we get

Thus the error limits are directly propor¬ tional to the standard deviation σ of the spectral noise and inversely proportional to the square root of the number of the data N. According to the definitions of inner product and norm defined by Eqs. 8 and 5, the remaining square root expression is not an explicit function of N. However, it depends on the shapes of the library spectra.

Let us now examine what happens when we dimin¬ ish the number of the data N by some factor 1/k. We see immediately the negative effect that the coefficient 1/ N is increased by factor k 12 . However, the standard deviation σ_ of the noise of the mixture spectrum has changed, too. This change is governed by the Parseval theorem

where n 1. and ns are the noise functions in the interfe- rogram and in the spectrum respectively. (These two random processes are a Fourier transform pair.) Because the noise n. is totally white, its "amplitude" is the

same everywhere. Thus when the length of the first integral is truncated to one k:th of its original va¬ lue, the value of the integral is diminished by the same factor k "1 . So the other integral has to change by same factor. Because we do not change the wavenumber range under investigation, the only possibility is that the "amplitude" of the noise n s , or its standard de- viation σ s , is diminished by factor k . This effect completely cancels the 1/ N dependence in Eq. 35. As we detected above, when reducing the reso¬ lution, the coefficient σ g /V " N on Eq. 35 remains con¬ stant. Thus the only possible source of changes in the error limits is the expression

As mentioned before, the definitions of inner product and norm mean that this expression depends solely on the shapes of the library spectra, when M is fixed. The number of the data itself is not important. All linear changes, where all the library spectra are multiplied by some constant coefficient C, change this square root expression by constant C "1 . Now the spectra are in prac¬ tice always computed from the corresponding interfero- grams by applying the Fast Fourier Transform (FFT) algorithm. A fundamental property of this algorithm is that the data interval in the spectrum is l/(2NΔκ), where Δκ is the sampling interval in the interferogram. So when the number of the data is diminished by factor 1/k, the data interval in the spectral domain is in¬ creased by factor k. As far as the data interval (=res- olution/1,21) stays smaller than the FWHH (full width at half height) of the spectral lines, there exists at least one data at every line and the shape of the spec- tral lines does not vary considerably. T i% means that there is only very little use from employing a better resolution than the width of the spectral lines. In the

interferogramn domain this means that we can safely truncate the interferogram provided that we do not cut off a considerable portion from the signal.

Let us now define the interferogram truncation function to be a boxcar function having the value of 1 between κ==pNΔκ and 0 elsewhere. Because we are only able to register a finite interval of the interfero¬ gram, the real, infinitely long interferogram is always multiplied by this function. In the spectral domain this means that the spectra are convolved with the Fourier transform w of the truncation function, or

e° = e "*w L , (37)

where w L (ι>) = 2NΔκsinc(2NΔκπυ), (38)

e° is the mixture spectrum or a library spectrum and e" is the spectrum which would be obtained by transforming the whole interferogram. The FWHH of this sine-function is approximately 1.21/(2NΔκ), and this is the entity we refer to as the resolution here. As far as N remains larger than (2ΔκxFWHH of the spectral lines) "1 , w L is narrower than the lines of e" and does not have a consi- derable influence on their shapes. If we still reduce the resolution after this point, the spectral lines suddenly begin to widen and their shapes become deter¬ mined mainly by w L instead of their real shapes. This means that the convolution of Eq. 37 then changes the spectra nonlinearly, so that its effect is not a mere multiplication of the square root expression by a cons¬ tant coefficient. If no apodizaton is applied, the lines begin to resemble sine-curves. (If apodization is performed, the interferogram noise data are no longer identically distributed and the error analysis is not valid.) Because of the widening, the lines begin to overlap, which makes it harder to discern them from

SUBSTITUTE SHEET

each other. This, in turn, means that the sum expressi¬ on in Eq. 35 begins to grow. The rate of this growth, however, depends on the amount M of the library spectra used in the analysis and on how near the lines are to each other. For example, if the lines were originally grouped into sets of overlapping lines, the growth rate would not be as fast as it would be were the lines originally situated at approximately equal intervals. Some crude results may be given anyhow. For instance, if we are using a set of at most 50 library spectra, the coefficient of growth of the square root expression normally lies between k 13 and k 1/2 depending on how the lines are situated. The square root expression also depends on the number of lines in the spectra. This dependence obeys the approximate law the he value of the square root is approximately inversely proportional to the square root of the average number of lines in one spectrum. Thus it can be considered as a constant coefficient not depending on N. According to what is said above, the best choice for the resolution would be the FWHH of the spectral lines. So the recorded interferogram should reach from -1/(2 x FWHH) to 1/(2 x FWHH). This, howe¬ ver, holds true only if we are unable to adjust the settings of the interferometer. If we can freely set all the parameters of the device, there exist two addi¬ tional advantages of the resolution reduction. We shall now consider them more closely.

As is well known, a non-point radiation source causes all spectral lines to widen. Typically the radi¬ ation source (opening) is round, and in that case every monochromatic spectral line is spread out into a boxcar line shape, as is shown in Fig. 3. The width of the box is then directly proportional to the area of the radia- tion source. This means that the spectrum e° ° in Eq. 37 is in fact the true spectrum e convolved with the box- cat function w arising from the non-zero area of the

light source. Thus Eq. 37 can be re-written as

e° = e*w Ω *w L . (39)

Because the width of w depends on the wavenumber υ Q of the spectral line under consideration, an exact treat¬ ment would require using different w for every line. An illustration of Eq. 39 is shown in Fig. 4. Because the FWHH of the convolved line is approximately the sum of the FWHHs of the components of the convolution, the distortions w D and w L have a considerable effect only, if their widths are larger than that of the natural widths of the spectral lines. Thus we can safely in¬ crease the signal by increasing the radius of the radi- ation source until the (maximal) width of the boxcar distortion equals the FWHH of the undistorted spectral lines. Correspondingly we can decrease the amount of computation by reducing resolution until the FWHH of the sine-distortion equals that of the spectral lines. (Which means the optimal truncation of the interfero¬ gram.) In case of gaseous specimen, however the natural width of the lines may be so small, that this situation can not be achieved. Anyhow, it is still sensible to make the distortions w π and w L of equal magnitude. Thus we can take this situation as our stating point. Let us now reduce the number of the data by factor k " . As we have stated before, this reduction widens the lines by widening w L by factor k thus increasing the square root expression of Eq. 35 by a factor of at most k 1/2, if we use some tens of library spectra. If M is of the order of a couple of hundreds, this factor can be of the order of k. Now, however, we can increase the area of the radiation source by factor k without any considera¬ ble further increase in the line widths. Because the area under the spectral lines must then grow by factor k due to the increase of the signal, the only possibi¬ lity is that the heights of the lines are increased by

the same factor k and the change of the spectra is approximately linear. Multiplication of the spectra with constant coefficient k reduces the square root expression by coefficient k "1 . This is more than enough to cancel the growth of the square root expression in the nonlinear interferogram truncation operation. In practice, however, there may arise difficulties in concentrating the enlargened image of the radiation source on the detector. The other additional advantage of the resolu¬ tion reduction is that we are now able to register k interferograms in the same amount of time we previously registered only one single. Because Fourier transfor¬ ming is a linear operation, co-adding these interfero- grams means that the corresponding spectra are co-ad¬ ded, too. The errorless spectra e remain the same in every measurement, which means that they become multip¬ lied by k in the summation. This means a simple linear change of the spectra, which in turn means that the square root expression becomes multiplied by k "1 . The noise of s, on the other hand is different every time, and it does not sum up linearly. From results 25 and 26 we see that the summed up noise has the distribution (N(0,kσ 2 ). Thus the standard deviation σ g of the noise is increased by factor k 1/2. The total effect is that the error limits become multiplied by factor k -1/2

When we finally gather together all the diffe¬ rent effects mentioned above, it becomes apparent that if we can freely change all the parameters of the in- terferometer, as small a resolution should be used as possible. The number of the data N should, however, be at least two or three times as large as the (maximal) number of the library spectra for the structure of the residual spectrum s-Kx to be examinable. If the library spectra are measured by using a different resolution than in measuring the mixture spectrum, the analysis can fail and large negative

coefficients may appear. A similar situation can arise if the line shapes in the library spectra and in the mixture spectrum differ due to nonlinearities or dif¬ ferent pressure broadenings. Then some improvement can be achieved by calculating the best non-negative solu¬ tion instead of the best solution. By non-negative solution we mean such coefficient vector x, which is the solution of the problem 4 subject to the condition that every component of x has to be non-negative. This procedure brings more information into the analysis, because we apply an a priori knowledge of the coeffi¬ cients. We are now going to derive an algorithm for finding the solution of problem 4 under the non-nega¬ tivity constraint. Let us denote d the residual norm

d(x) = s-Kx

Because a norm is always a non-negative quantity, the norm ||d|| has exactly the same minima as its square ||d|| 2 so that instead of the norm we can minimize its square. Now ||d(x)|| 2 is a convex function of x. This means that for every x χ , X 2 and λ, 0<λ<l,

|| d(Ax x + (1 - λ)x 2 ) || 2 < λ || d( Xl ) +(1 - λ) || d(x 2 ) || 2 .

This can be seen by using the triangle inequality and the fact that the geometrical mean is always less or equal to the arithmetic mean. The convexity implies that the square norm has only one minimum point, which makes the minimization much simpler.

Especially when other coefficients are held fixed and only one coefficient κ, is varied, ||d|| is a convex function of one variable. Thus at the optimum point x there exist two possibilities. Either

d || d(x opt ) || 2 _ Q dxj

S

(where all the components except the j : th are fixed to xopt)' or,' if the zero of he derivative does not lie in the allowed area κ 3>0,

5 κ o P t, j = °'

which means that κopt, _. is on the border between the allowed and the forbidden areas. This can be proved as 10 follows:

1) If the zero of θ||d(x ot )|| /θκ lies in the allowed area κ>0, κ __, clearly has to equal this zero.

2) If the zero of d||d(x o t ) || 2 /θκ lies in the forbidden area κ <0, the derivative is positive, when κ > 0,

15 because of the convexity of ||d(x)|| 2 . Thus, if κ o > 0, diminishing of κ o would diminish the value of ||d|| 2 without exiting the allowed area. Thus the only possi¬ bility 1 is that

We are now able to state the following condi- 20 tion for optimality:

At the unique solution point x of the minimization problem 4 with non-negative component" x either of the following conditions holds for eac-, .ponent κ.,:

25 0 and d » «g: > o) . (

This is illustrated in Fig. 5. The so-called Kuhn-Tucker criteria (see some textbook oi optimizati¬ on, e.g. [3]), which in case of a convex object functi- 30 on are a necessary and a sufficient condition for opti¬ mality, would give quite the sane condition.

The partial derivatives are - rγ simple to calculate, and we obtain

3.5 V χ (||d(x)|| 2 ) = 2(Ax-y), (41)

where A is the inner product matrix defined in Eq. 10

and y is given by Eq. 11. The individual partial deri¬ vatives are the components of this gradient, or

5(11d pII 2 )-=(V x ||d(x)|| 2 )-=2 (T A jiXi - yή (42) dx

If we include some background spectra in the used set of library spectra, it is sensible not of impose the non-negativity constraint on their coeffi- cients. The background spectra consist of all the com ¬ pounds, which are present also when measuring the back¬ ground (such as water and carbon dioxide), plus a con¬ stant function and simple cosine curves, which can be generated by some erroneous data in the interferogram. When this is taken into account, too, the minimum point can be found e.g. by the following algorithm: 1) Choose some starting point x=x 0/ for instance x =(0,0, ... ,0) τ or x = A _1 y Insert j=l.

2) Calculate

Xj _ = y j - ∑i≠j -Aji X i

~-jj which according to Eq. 42 is the zero point of

|2 d d(x /dx . . If κ.<0 and the coefficient x_ has the non-

11 " 3 3 j negativity constraint, re-insert κ =0. Insert j=j+l.

3) If j≤M (where M is the amount of library spectra used in the analysis), return to step 2) . Otherwise continue form step 4).

4) Calculate

G = iv x (||d(x)|| 2 ) = Ax-y.

If for every component G of G either | G | <e. or (G > 0 and κ < e 2 ), where e χ and e 2 are suitable small real numbers, stop. Otherwise insert j = 1 and return to step 2) .

References

1. C. L. Lawson & R. J. Hanson, Solving least squares problems (Prentice-Hall, Ne Jersey, 1974), p. 36.

2. B. Noble <§c J. W. Daniel, Applied linear algebra (Prentice-Hall, New Jersey, 1977 2nd ed.

3. M. S. Bazaraa & C. M. Shetty, Nonlinear programming - theory and algorithms (Joh Wiley & Sons, New York, 1979), p. 137.