Efficient and robust method for simultaneous reconstruction of the temperature distribution and radiative properties in absorbing, emitting, and scattering media

Efficient and robust method for simultaneous reconstruction of the temperature distribution and radiative properties in absorbing, emitting, and scattering media

Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57 Contents lists available at ScienceDirect Journal of Quantitative Spectro...

2MB Sizes 1 Downloads 44 Views

Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

Contents lists available at ScienceDirect

Journal of Quantitative Spectroscopy & Radiative Transfer journal homepage: www.elsevier.com/locate/jqsrt

Efficient and robust method for simultaneous reconstruction of the temperature distribution and radiative properties in absorbing, emitting, and scattering media Chun-Yang Niu, Hong Qi n, Xing Huang, Li-Ming Ruan, He-Ping Tan School of Energy Science and Engineering, Harbin Institute of Technology, Harbin 150001, PR China

a r t i c l e in f o

abstract

Article history: Received 18 January 2016 Received in revised form 21 June 2016 Accepted 22 June 2016 Available online 29 June 2016

A rapid computational method called generalized sourced multi-flux method (GSMFM) was developed to simulate outgoing radiative intensities in arbitrary directions at the boundary surfaces of absorbing, emitting, and scattering media which were served as input for the inverse analysis. A hybrid least-square QR decomposition–stochastic particle swarm optimization (LSQR–SPSO) algorithm based on the forward GSMFM solution was developed to simultaneously reconstruct multi-dimensional temperature distribution and absorption and scattering coefficients of the cylindrical participating media. The retrieval results for axisymmetric temperature distribution and non-axisymmetric temperature distribution indicated that the temperature distribution and scattering and absorption coefficients could be retrieved accurately using the LSQR–SPSO algorithm even with noisy data. Moreover, the influences of extinction coefficient and scattering albedo on the accuracy of the estimation were investigated, and the results suggested that the reconstruction accuracy decreased with the increase of extinction coefficient and the scattering albedo. Finally, a non-contact measurement platform of flame temperature field based on the light field imaging was set up to validate the reconstruction model experimentally. & 2016 Elsevier Ltd. All rights reserved.

Keywords: Generalized sourced multi-flux method Inverse radiation problem LSQR–SPSO algorithm Retrieval temperature field

1. Introduction Participating media are widely presented in many engineering fields, such as aerospace engineering, energy and power systems, ferrous metallurgy, chemical engineering, and information communications [1–6]. In these fields of practical interest, the temperature distribution and radiative properties of participating media all play an important and indispensable role in the analysis of heat transfer or information transfer [1–6]. Past years have witnessed sustained efforts aimed at measuring the temperature distribution and the radiative properties of participating media [5–8]. Various measuring techniques have n

Corresponding author. Tel./fax: þ86 0451 8612638. E-mail address: [email protected] (H. Qi).

http://dx.doi.org/10.1016/j.jqsrt.2016.06.032 0022-4073/& 2016 Elsevier Ltd. All rights reserved.

been developed to obtain these parameters. In general, the temperature measurement methods of participating media can be classified into two categories, namely, the contact measurement method and non-contract (or nonintrusive) measurement method. Contact measurement methods have the disadvantages of perturbing flow-field by physical probing, causing measurement errors by physical and chemical reactions that may occur within a sampling manifold [9]. So these method cannot be used in many high-accuracy required measurement situations, especially for measuring the three-dimensional (3D) temperature field. However, the non-contact optical measurement technology based on the radiation image processing is non-invasive to cause minimal disturbance of the media while the media is being probed, and can be used for in-situ measurement with high precision. Especially, the information of the 3D temperature distribution

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

Nomenclature Fobj I Ib In In1

the the the the the

In2

the measured radiative intensity in direction

0

In2 L M n SNR R s S T T 0 T X0,Y0,Z0

ψ1

objective function radiative intensity [W/(m2 sr)] blackbody radiation intensity [W/(m2 sr)] exit radiative intensity vectors measured radiative intensity in direction

ψ2

the assumed radiative intensity in direction ψ2 the length of the medium [m] the number of detection line the grid number along the detection line the signal-to-noise ratio the radius of the medium [m] the location the generalized sourced the temperature [K] the temperature field [K] the assumed temperature field [K] the start point of the detecting line

also can be obtained through this method. Thus, the noncontact optical measurement technology has been recognized as a potentially efficient mean of accurately monitoring the temperature of flame in combustion systems [5–9]. The estimation accuracy of temperature distribution is directly related to the radiative properties of participating media, and more accurate data with radiative properties would be of significant help to reconstruct the temperature distribution accurately. The problems of predicting the temperature distributions and radiative properties from the radiative energy images can be merged into an inverse radiation problem. The investigation of retrieval temperature distribution problems for participating media of complex radiative heat transfer is far from complete. Nowadays, problems in retrieval temperature distribution have drawn considerable attention in the field of combustion and represent a highly promising research direction and potential practical technology. The core problem of non-contact optical measurement technology based on radiation image processing is to develop an efficient and accurate forward solving method for simulating the outgoing radiative intensities in arbitrary directions at boundary surfaces. An efficient, accurate, and robust inverse problem-solving mode is urgently needed to retrieve the temperature distribution and radiative properties of participating media. The general challenges for this task are concerned with mathematical ambiguity. Most of the time, some simplifying assumptions (e.g., the assumption of pure absorbing participating medium by ignoring the scattering effect) are posited to reduce the difficulties; however, the simplifications must be reasonable for realistically modeling the physical processes [5]. During the last two decades, considerable researches have been conducted to determine the temperature profile or source term and radiative properties in

45

Greek symbols

ε Φ κa κe κs ω θ ϕ Ω ς

the the the the the the the the the the

error scattering phase function absorption coefficient [m  1] extinction coefficient [m  1] scattering coefficient [m  1] scattering albedo polar angle of the detection direction azimuthal angle of the detection direction solid angle [sr] normal distribution random variable

Subscripts est exa mea mean rel

the the the the the

estimated value exact value measurement value mean value relative value

participating media from various types of radiation measurements [9–19]. A comprehensive review of inverse radiative problems in participating media has been given by McCormick [20]. Most studies have assumed the radiative properties to be independent of wavelength because of the lack of reliable and sufficient experimental data. Most of the analytical or numerical methods dealing with the inverse radiation problems of estimating the temperature fields are one- or two-dimensional (1D or 2D). Only a few works estimate the unknown temperature sources in 3D problems because of the complexity of geometry or time-consuming computation. For instance, Liu and Tan et al. [1] developed a hybrid network searchconjugate gradient algorithm for simultaneously estimating temperature distribution and boundary emissivities for an emitting, non-scattering, gray, and plane–parallel medium from the knowledge of the exit radiation intensities at boundary surfaces. Liu et al. [21] also investigated the simultaneous identification of temperature profile and absorption coefficient in 1D semitransparent media. Zhou et al. [11] have demonstrated the potential to simultaneously determine the temperature distribution and radiative properties for 1D media using the Newton-type least-square iteration algorithm. Lou et al. [22] proposed a decoupled method for simultaneously estimating temperature distribution and the radiative properties for 1D optical thin axisymmetric flame or non-scattering participating media. Liu et al. [23] developed an inverse radiation analysis to simultaneously estimate temperature field and uniform radiative properties in participating medium using least-square QR decomposition (LSQR) method, combined with the backward Monte Carlo (BMC) method, and used charge-coupled device (CCD) cameras to simultaneously estimate 3D soot temperature and volume fraction fields in optically thin flames [24,25]. Huang et al.

46

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

[26] established a new iterative reconstruction model for simultaneously estimating 3D soot temperature and volume fraction distributions in unsteady sooty flames. Snelling et al. [27] developed a multi-wavelength flame emission technique based on the line-of-sight (LOS) method for the simultaneous estimation of soot temperature and volume fraction in axisymmetric laminar diffusion flames. More recently, Daun et. al. [28,29] employed infrared laser light absorption measurement to retrieve the chemical species concentration distribution within in combustion devices. The Tikhonov regularization based on a Laplacian smoothing matrix is adopted to solve the limited data chemical species tomography (CST) problem and the level set method is used to solve the line-of-sightattenuation chemical species tomography [28]. Furthermore, various ways that prior information can be used to enhance reconstruction accuracy of CST experiments on turbulent flows are reviewed more recently [29]. As a whole, most of the works on retrieving temperature field in participating media must know (or assume) the radiative properties of the media beforehand. Moreover, some of the works on the simultaneous estimation of 3D temperature distribution and radiative properties in participating media also need to assume the media is pure absorbing media. In fact, accurate determination of the radiative properties of dispersed media should be considered today as an unsolved problem that requires deep investigations. Further progress in reconstructing 3D temperature distribution is necessary to simultaneously estimate the radiative properties of absorbing and scattering media. Nevertheless, to the best of the authors' knowledge, only a few studies have focused on the simultaneous estimation of 3D temperature distribution and radiative properties in absorbing, emitting, and scattering media based on solving the direct problem of radiative transfer. To date, the Monte Carlo method (MCM) is used as the direct problem-solving method for simulating outgoing radiative intensities in arbitrary directions at the boundary surfaces of the absorbing, emitting, and scattering media in most works. Although the MCM has high calculation accuracy, it requires heavy computation efforts, and its calculation speed is extremely slow to affect the efficiency of the overall temperature reconstruction algorithm. For practical purposes, a novel direct problemsolving method and corresponding inverse algorithm with high precision and efficiency is urgently required to simultaneously reconstruct 3D temperature distribution and radiative properties in absorbing, emitting, and scattering media. In this study, a rapid and accurate computation method called generalized sourced multi-flux method (GSMFM) is developed to simulate outgoing radiative intensities in arbitrary directions at the boundary surfaces of the absorbing, emitting, and scattering media. A hybrid LSQR– stochastic particle swarm optimization (LSQR–SPSO) algorithm is also developed as an inverse method for simultaneously reconstructing the multi-dimensional temperature distribution, absorption coefficient, and scattering coefficient theoretically and experimentally. The influences of measurement errors and radiative properties on the reconstruction accuracy of temperature field and radiative

properties are also investigated. The remainder of the paper is organized as follows. The solving model for direct problem is examined in Section 2. The inversion principle of the LSQR–SPSO algorithm is described in Section 3. The retrieval results of temperature field and radiative properties in participating media using the LSQR–SPSO algorithm and the effects of the measurement errors on the accuracy of the inverse analysis are presented in Section 4. An experimental platform is developed to validate the reconstruction model in Section 5. The main conclusions and perspectives are provided in Section 6.

2. Direct problem 2.1. Physical model The direct problem solver is particularly critical for an optimization approach. For the simultaneous reconstruction of the temperature distribution and radiative properties in absorbing, emitting, and scattering media, the first and most essential stage is to develop an efficient and accurate forward solving method. Considering the cylindrical participating media, as shown in Fig. 1, all of the boundary surfaces are supposed to be transparent. The radiation intensity in arbitrary direction can be solved using the GSMFM. It is a rapid computation method based on the finite volume method (FVM) for the exact solution of the arbitrary directional radiative intensity of the participating media. It takes into account computational efficiency and accuracy comprehensively, thereby achieving the computational accuracy close to the MCM and obtaining the computational efficiency and applicability close to the FVM. For the participating media, the radiative transfer equation (RTE) takes the following form: dIðs; sÞ ¼  κ e I ðs; sÞ þ Sðs; sÞ ds Z κs i I ðs; si ÞΦðsi ; sÞdΩ Sðs; sÞ ¼ κ a I b ðsÞ þ 4π 4π

ð1Þ ð2Þ

where I(s, s) represents the radiative intensity at position s and direction s; Ib(s) denotes the blackbody radiation

Fig. 1. The physical model of cylindrical participating media

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

intensity; κa, κs and κe are the absorption coefficient, the scattering coefficient, and extinction coefficient of the media, respectively; Φ(si, s) denotes the scattering phase function; and S(s, s) represents the generalized source term. According to the FVM model, the radiation energy conservation equation can be expressed as M X j

Z Aj Ij ðsi Þ

 ΔΩ

i

 i i si U nj dΩ ¼ ½  κ e Iðs; si Þ þ Sðs; si ÞΔV ΔΩ

i ¼ 1; ⋯; N s

ð3Þ where Aj is the surface area of control volumes; M is the surface number of control volumes; nj is the normal vector of the outside surface of control volumes; and Ns is the number of discrete directions. To simulate the radiative intensity of detection direction sk, the new direction sk must be added into the total discrete directions (Ns þ1). For the solid angle in the simulated direction sk is extremely small (ΔΩk E0), Eq. (3) can be approximated as M X j

Z Aj I j ðsk Þ

ΔΩ

k

  k k sk Unj dΩ ¼ ½  κ e I ðs; sk Þ þ Sðs; sk ÞΔV ΔΩ  0

47

The integral of Eq. (6) in a micro-control volume is obtained as follows:   I i þ 1 BA ¼ I i BA expð  κ e Δsi Þ þ 1  expð  κ e Δsi Þ SBA ð7Þ i SBA i ¼ ð1  ωÞI bi þ

N ωX

4π j ¼ 1

  s j I i j Φ sj ; sBA ΔΩ

ð8Þ

The generalized source Si can be obtained using FVM. The boundary surfaces of the participating media are supposed to be transparent. Therefore, if the environmental radiation is ignored, the exit radiative intensity in any direction r can be obtained by piecewise integrating Eq. (7) from the boundary surfaces of the media to the detector. Hence, the equation can be written as follows: 2 0 1 1 n X   nX BA  BA 4exp@  I ¼ S 1  exp  κ e Δsi þ κ e Δs j A n

n

0  exp@ 

n X

13

i¼1

κ e Δsj A5SBA i

j ¼ iþ1

ð9Þ

j¼i

The detailed solution procedure of the GSMFM is available in Ref. [30] and not repeated here.

ð4Þ The linear equations number of the original FVM model will not increase, which signifies that new direction sk has no influence on the original FVM calculation. Hence, Eq. (2) can be written as Sðs; si Þ ¼ κ a I b ðsÞ þ

Ns     κs X j I s; sj Φ sj ; si ΔΩ i ¼ 1; ⋯; Ns þ 1 4π j ¼ 1

ð5Þ Therefore, the basic idea of GSMFM is to calculate the directional radiation intensity in the Ns discrete directions by solving the RTE using FVM firstly and then calculating the generalized source S(s,s). As shown in Fig. 2, RTE can be established along the arbitrary direction of BA to calculate the directional radiation intensity in the direction of BA, which can be expressed as dI ðs; sBA Þ ¼  κ e I ðs; sBA Þ þ SBA s ds

ð6Þ

2.2. Validation of the GSMFM model Given the cylindrical participating media model, as shown in Fig. 1, the radius is R¼0.5 m and the axial length is L¼20 m. The cylindrical media is divided into 22880 grid, that is, Nϕ  Nr  Nz ¼20  22  52. The solid angle is divided as Nθ  Nϕ ¼10  9. The temperature of the media is set as T¼ 2000 K. The scattering phase function is assumed as Φ ¼1, the extinction coefficient is κe ¼ 2.0 m  1, and the scattering albedo is ω ¼0.2, 0.5, 0.8. Given a set of detecting line, the coordinate of start point of the detecting line (x0, y0, z0) is set as ( 0.5,  0.5, 10) [m]. The azimuthal angle is ϕ ¼45°, and the polar angle θ is in the range of [10–170°]. The outgoing radiative intensities of the set of detecting lines are calculated using GSMFM and BMC, respectively. The simulated results are illustrated in Fig. 3. As shown in Fig. 3, the results calculated by GSMFM are consistent with those calculated by BMC. The maximum relative error between GSMFM and BMC is about 2.3%. Meanwhile, the computing time of the GSMFM and BMC algorithm under the same condition are shown in Table 1. The computation speed of the GSMFM is faster than that of the BMC. The computational time of the GSMFM algorithm is only about 1% of the BMC algorithm, which proves that GSMFM has high computing efficiency. Thus, GSMFM can be regarded as a reasonable trade-off between efficiency and accuracy, which renders the suitability of GSMFM as a forward calculation method for solving the inverse problem.

3. Inverse principle

Fig. 2. Computational path of the directional radiative intensity

In the inverse problem, the radiative intensities received by the detectors are the measurement signals, and the temperature field and the radiative properties are

48

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

Fig. 3. The outgoing radiative intensities in different directions Table 1 The computing time of the GSMFM and BMC algorithm. Method

GSMFM

BMC

Time/s

51

4850

the unknown parameters. According to Eqs. (8) and (9), the temperature field and radiative properties can be retrieved from the knowledge of the outgoing radiative intensities. In fact, the combined reconstruction of the temperature field and radiative properties is a typical illposed inverse problem. For the generalized source and radiative properties are coupled together when the radiative properties is unknown from Eqs. (1) and (2), it makes the integral radiative transfer equations along the detection lines (see Eq. (9)) to be nonlinear. Then the reconstruction problem of temperature distribution and radiative properties becomes a nonlinear optimization problem. The regularization methods (such as the LSQR algorithm) cannot solve this problem independently. The intelligent optimization algorithms (such as SPSO) and the gradient-based method (such as CG) can solve the nonlinear optimization problem, but they are powerless for solving the reconstruction problem of large-scale multidimensional and multiple parameter field independently. Therefore, a hybrid LSQR–SPSO method is developed based on the least-square QR decomposition algorithm and the SPSO algorithm in present research. The SPSO algorithm is used to retrieve the radiative properties with known retrieval temperature distribution by LSQR, and the LSQR algorithm is used to solve the temperature distribution with known retrieval radiative properties by SPSO. On this basis, the multi-dimensional temperature distribution and radiative properties can be reconstructed simultaneously by the hybrid LSQR–SPSO algorithm. It is worth noting that the drawbacks to using the metaheuristics like the LSQR– SPSO algorithm are: (1) searching for multiple local minima that don't exist is inefficient and time consuming; (2) there is no guarantee that the final point is stationary (i.e. the gradient is zero) unless there is a nonlinear

programming ‘finishing’ step such as conjugate gradient method to obtain the minimum. These drawbacks needs to be improved in the further study. The LSQR algorithm was first developed by Paige and Saunders [31] in 1982. Its basic idea is to turn the arbitrary coefficient matrix equations into equations whose coefficient matrix is a square matrix at first, and to subsequently solve the least-square solution of the equations using the Lanczos method. The LSQR algorithm has quick convergence, high computing precision, and good numerical stability, and is suitable for solving an ill-posed problem. It is particularly suitable for solving linear equations with a large-scale sparse coefficient matrix. The PSO algorithm was first developed by Eberhart and Kennedy [32] in 1995. It is an optimization algorithm based on swarm intelligence with random search. Its basic idea is derived from the foraging behavior of birds. Group cooperation and competition among the particles produce swarm intelligence to guide the optimization search. The PSO algorithm, as an efficient and simple parallel search algorithm, has been widely applied in all types of optimization problems in recent years [33–36]. Qi and Ruan et al. [10,33] introduced the PSO algorithm into the inverse problem of radiation source and radiative properties and reported the efficiency and robustness of the PSO algorithm for inverse radiation problems with the advantages of being independent from the initial value and eradicating the need to solve the objective function derivative. The SPSO algorithm was first developed by Zeng et al. [37]. Compared with the standard PSO, the SPSO eliminates the historical velocity term, which makes the particle lose velocity memory to decrease the global searching capability. However, SPSO algorithm guarantees that, in every generation, one particle stops evolving so that best position could be located. The basic concept of the SPSO algorithm is stopping the particle from evolving to improve the global searching capability. The SPSO algorithm could guarantee global convergence, whereas the PSO could not. A hybrid LSQR–SPSO algorithm is developed by combining the stochastic searching capability of the SPSO algorithm and the ill-posed inverse problems solving capacity of the LSQR algorithm. The objective function is defined as F obj ðT; κ a ; κ s Þ ¼

‖I0n2  In2 ‖ kIn2 k

ð10Þ

where In2 is the measured radiative intensity in direction 0 ψ2, which is calculated by the direct problem, and In2 is the assumed radiative intensity of direction ψ2 with the 0 assumed temperature field T . The value of objective function under different absorption coefficient and scattering coefficient is calculated. The results are shown in Fig. 4. The absorption coefficient and scattering coefficient are set as κa ¼1.0 m  1 and κs ¼3.0 m  1, respectively. From Fig. 4, it can be seen that the value of objective function will be minimum only when the absorption coefficient and scattering coefficient are close to the exact value, which proves the objective function appears convex.

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

49

Fig. 4. The value of objective function under different absorption coefficient and scattering coefficient

Moreover, the implementation of the LSQR–SPSO approach for solving the inverse radiative problem can be performed according to the following procedure: Step 1. The measured radiative intensities In1 and In2 in two different directions ψ1 and ψ2 calculated by GSMFM are selected as input data for simultaneously retrieving the temperature distributions, absorption coefficient, and scattering coefficient. Step 2. Systemic parameters of SPSO are input and the particle swarm is initialized, where the particle position represents the absorption coefficient and scattering coefficient of the media. 0 Step 3. The assumed generalized source S is calculated by the LSQR algorithm from the measured radiative intensity In1 of direction ψ1 and the particle position (i.e. absorption coefficient and scattering coefficient). 0 Step 4. The assumed temperature field T is calculated 0 by FVM according to the assumed generalized source S . 0 Step 5. The assumed radiative intensity In2 of direction ψ2 is calculated by GSMFM from the assumed temperature 0 field T and the particle position (i.e. absorption coefficient and scattering coefficient). The fitness value is subsequently calculated by Eq. (10).

Fig. 5. The flowchart of the hybrid LSQR–SPSO algorithm

Step 6. For each particle, the experienced fitness values are compared to determine the current individual best position. The current individual best position is consequently selected as the current global best position.

50

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

Fig. 6. The simplified model of simulate light field imaging

Step 7. For each particle, the new velocity is calculated and the new position is updated. Step 8. The stop criterion is checked. If the fitness value is less than the setting tolerance or the iteration number exceeds the maximum number, then proceed to Step 9; otherwise, loop to Step 3. Step 9. The retrieval results of temperature field Test are calculated using the LSQR algorithm from the measured radiative intensity In1 of direction ψ1 and the global best position (the retrieval absorption coefficient κa, est and the retrieval scattering coefficient κs, est). Step 10. The optimal solution of the inverse radiative problem is generated. The detailed flowchart of the simultaneous retrieval of the multi-dimensional temperature distributions and radiative properties using the hybrid LSQR–SPSO algorithm is presented in Fig. 5.

4. Results and discussion Given the simplified model of simulated light field imaging [19], as shown in Fig. 6, the radiation transmitted along straight lines in the actual imaging equipment and the extinction of the radiation in the imaging equipment is ignored, only the radiative intensities on the pupil plane of imaging equipment are calculated in the numerical simulation. Therefore the four-dimensional (4D) light field (include 2D position information and 2D direction information) is established on the pupil plane of the imaging equipment, which is divided into 10  10 macro pixel grids (equivalent to the amount of micro lens), and each macro pixel grid is divided into 20  20 detection directions, the total number of the detection lines is M ¼40,000. Moreover, the radius is R¼0.5 m, and the axial length is L ¼20 m. The cylindrical media is divided into 300 grids, that is, Nϕ  Nr  Nz ¼ 3  10  10. To demonstrate the effects of measurement noise on the predicted terms, the measured radiative intensity with random noise is obtained by adding normal distribution

noise to the exact one, as shown below:   I mea ¼ 1 þ φς I exa

ð11Þ

where I represents the radiative intensity exiting the boundaries, and ϛ is a normal distribution random variable within a range of  2.576 to 2.576 with the probability of 99%. ϕ is the mean square deviation of the noise, which can be determined as SNR ¼ 20 logðφÞ

ð12Þ

where SNR represents the signal-to-noise ratio. For the sake of comparison, the relative error is defined as   ϵrel;i ¼ T est;i  T exa;i =T exa;i ð13Þ where Test,i represents the retrieval temperature, and Texa,i denotes the exact temperature. The mean relative errors of retrieval results of the temperature field are defined as

ϵrel;mean ¼

N 1X ϵ N i ¼ 1 rel;i

ð14Þ

where N is the grid number of the cylindrical media. Case 1. Axisymmetric temperature distribution. The exact axisymmetric temperature field of the media is assumed to be the function of z and r as given by T ðz; r Þ ¼ 

  2600 ðz  20:0Þ  0:5 3:0r 2 þ 2:0r 3 þ 800 17

ð15Þ

where z represents the axial coordinate of the cylinder, and r is the radial coordinate of the cylinder. The temperature distribution is shown in Fig. 7. In this case, the exact values of the absorption and scattering coefficients are κa ¼1.0 m  1 and κs ¼ 3.0 m  1, respectively. The searching ranges of κa and κs are all set as [0.1, 4.0] m  1. The SNR of the measurement is selected as 40, 60, and 70 dB. The 2D temperature distribution, absorption coefficient, and scattering coefficient are reconstructed simultaneously using the hybrid LSQR–SPSO algorithm. The radiative intensities received by the detectors at two detection directions of θ ¼90°, ϕ ¼30°

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

51

r

r

0.5

0.0

0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20

-0.5 0

2

4

6

8

10 z

12

14

16

18

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

No noise Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

r

Fig. 7. The exact axisymmetric temperature distribution of media

Table 2 The retrieval results of radiative properties in Case 1.

No noise 70 dB 60 dB 40 dB

1

κa [m ] κs [m  1] κa [m  1] κs [m  1] κa [m  1] κs [m  1] κa [m  1] κs [m  1]

True value

Retrieval value

εrel %

1.0 3.0 1.0 3.0 1.0 3.0 1.0 3.0

0.9999 3.0001 0.9997 3.0003 1.0041 2.9959 0.9773 3.0230

0.01 0.00 0.03 0.01 0.41 0.14 2.27 0.77

and θ ¼60°, ϕ ¼ 45° are selected as input measurement data. The retrieval results of the absorption and scattering coefficients are listed in Table 2. The retrieval temperature field with different SNR values is shown in Fig. 8. The relative reconstruction error distributions of temperature field with different SNR values are presented in Fig. 9. The mean reconstruction error of the temperature field and the reconstruction error of radiative properties under different SNR values are shown in Fig. 10. As shown in Table 2 and Fig. 10, the absorption and scattering coefficients could be retrieved accurately without measurement noise. The reconstruction error increases with the decrease of the SNR value. The maximum reconstruction error of the absorption coefficient is only 2.27%, and the maximum reconstruction error of the scattering coefficient is merely 0.77% even with SNR of 40 dB. Therefore, the absorption and scattering coefficients could be retrieved reasonably using the hybrid LSQR–SPSO algorithm even with noisy data. As shown in Figs. 8–10, the retrieval temperature field is accurate in case of no measurement noise, and the relative error is in the range of 10  5 to 10  6. The reconstruction error increases with the decrease of the SNR value, but the maximum reconstruction error is merely 1.88% even with SNR¼40 dB, and the magnitude of mean reconstruction error is only in the order of 10  3. Thus, the 2D axisymmetric temperature field can be retrieved

0

5

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

SNR=70dB

r

Inverse parameter

0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

SNR=60dB

r

SNR

0.5 0.4 0.3 0.2 0.1 0

0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

SNR=40dB Fig. 8. The retrieval results of temperature field in Case 1

accurately using the hybrid LSQR–SPSO algorithm, even with measurement errors, which proves that the hybrid LSQR–SPSO algorithm is accurate and robust. Case 2. Non-axisymmetric temperature distribution. Not all temperature fields of media are axisymmetric in practical engineering applications. A large part of the temperature field is typically irregular. Therefore, the reconstruction of the non-axisymmetric temperature field

     2 2 T ðx; y; zÞ ¼ 1000 U exp  40 ð12x þ 7:5Þ=9 1:1 25 ð12y þ8:5Þ=9  0:8    2  2 þ0:8exp 25 ð12x þ 7:5Þ=9  0:8  35 ð12y þ 8:5Þ=9  1:2 þ 60 U ð20  zÞ þ 800

ð16Þ

52

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

No noise

0

50

100

150

200

250

300

Elements

NSR=70dB

0

50

100

150

200

250

300

Fig. 10. The mean reconstruction error distribution of temperature field and radiative properties in Case 1

Elements

Tr/K

0.5

3000 2950 2900 2850 2800 2750 2700 2650 2600 2550 2500 2450 2400 2350 2300 2250 2200 2150 2100 2050 2000

0.4 0.3 0.2

NSR=60dB 50

100

150

200

250

0.1

300

y

0

Elements

0 -0.1 -0.2 -0.3 -0.4 -0.5 -0.5 -0.4 -0.3 -0.2 -0.1 0 x

NSR=40dB 0

50

100

150

200

250

0.1 0.2 0.3 0.4 0.5

300

Fig. 9. The reconstruction error distribution of temperature field in Case 1

and radiative properties are required to be examined. The non-axisymmetric bimodal temperature distribution is given by where x, y, z represent the coordinates of the cylindrical media. The temperature distributions in the sections of z ¼0 and ϕ ¼300° are shown in Fig. 11(a) and (b), respectively. In this case, the exact values of the absorption and scattering coefficients are also set as κa ¼1.0 m  1 and κs ¼3.0 m  1. The searching ranges of κa and κs are all set as [0.1, 4.0] m  1. The SNR of the measurement is assumed to be 40, 60, and 70 dB. The temperature distribution, absorption coefficient, and scattering coefficient are reconstructed simultaneously using the hybrid LSQR–SPSO algorithm. The radiative intensities received by the detectors at two detection directions of θ ¼90°, ϕ ¼30°

r

Elements

0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300 Fig. 11. The temperature distribution in the section of z ¼0 (a) and the section of φ ¼ 3001 (b)

and θ ¼60°, ϕ ¼45° are selected as input. The retrieval results of the absorption and scattering coefficients under different SNR values are listed in Table 3. The retrieval results of the temperature field with different SNR values are shown in Fig. 12. The relative reconstruction error distributions of temperature for all elements under different SNR values are shown in Fig. 13. The mean

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

53

Table 3 The retrieval results of radiative properties in Case 2. SNR

Inverse parameter

No noise 70 dB 60 dB

r

40 dB

0.5 0.4 0.3 0.2 0.1 0

0

κa κs κa κs κa κs κa κs

1

[m ] [m  1] [m  1] [m  1] [m  1] [m  1] [m  1] [m  1]

5

True value

Retrieval value

εrel %

1.0 3.0 1.0 3.0 1.0 3.0 1.0 3.0

0.9992 3.0008 0.9983 3.0017 1.0057 2.9947 1.0976 2.9024

0.08 0.03 0.17 0.06 0.57 0.18 9.76 3.25

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

r

No noise 0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

r

SNR=70dB 0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20 Fig. 13. The reconstruction error distribution of temperature field in Case 2

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

r

SNR=60dB 0.5 0.4 0.3 0.2 0.1 0

0

5

10 z

15

20

Tr/K: 900 1040 1180 1320 1460 1600 1740 1880 2020 2160 2300

SNR=40dB Fig. 12. The retrieval results of temperature field in Case 2

reconstruction errors of the temperature field and the reconstruction errors of radiative properties under different SNR values are shown in Fig. 14.

Fig. 14. The mean reconstruction error distribution of temperature field and radiative properties in Case 2

54

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

10-1

Reconstruction error

Reconstruction error

Mean of temperature Absorption coefficient Scattering coefficient

10-2

10-3

10

-4

10

-5

10-6

10

-1

10

-2

10

-3

10

-4

-1

10

0.25

0.5 0.75 Scattering albedo ω

Mean of temperature Absorption coefficient Scattering coefficient

10-5

κe =4.0 m

0

ω =0.667

1

-6

1.5

3 Extinctioncoefficient κe / m-1

4.5

Fig. 15. The effects of the scattering albedo on the reconstruction accuracy of temperature field and radiative properties

Fig. 16. The effects of the extinction coefficient on the reconstruction accuracy of temperature field and radiative properties

As shown in Table 3 and Fig. 14, the absorption and scattering coefficients could also be retrieved accurately without measurement errors. When measurement error exists, the reconstruction error increases with the decrease in the SNR value. The results indicate that the absorption and scattering coefficients could be retrieved reasonably by the hybrid LSQR–SPSO algorithm even with noisy data. As shown in Figs. 12–14, the retrieval temperature field is accurate for the case without measurement errors, and the relative error is in the range of 10  3 to 10  6. The reconstruction error increases with the decrease of the SNR value, but the maximum reconstruction error is only 2.1% even with SNR¼40 dB, and the magnitude of mean reconstruction error is merely in the order of 10  2. Thus, the non-axisymmetric temperature field can also be retrieved accurately using the hybrid LSQR–SPSO algorithm, even with measurement errors. In addition, compared with the radiative properties, the temperature field is easier to be retrieved and exhibits higher accuracy because the radiative intensity is insensitive to the absorption and scattering coefficients, which may lead to the multi-value characteristic of the retrieval radiative properties. In other words, different multi-retrieval results of the absorption and scattering coefficients all can make the objective function equal to an extremely small value. The effects of the extinction coefficient and scattering albedo on the reconstruction accuracy of temperature field as well as the absorption and scattering coefficients are also investigated in the case of non-axisymmetric bimodal temperature distribution. First, the exact extinction coefficient of the medium is assumed to be 4.0 m  1, and the exact scattering albedo is set as 0.05, 0.25, 0.5, 0.75 and 0.95 to retrieve the temperature distribution and radiative properties using the hybrid LSQR–SPSO algorithm. The mean reconstruction error of the temperature field and the

reconstruction error of radiative properties are shown in Fig. 15. Consequently, the exact scattering albedo is assumed to be 0.667, and the exact extinction coefficient of the media is set as 1.5, 3.0, and 4.5 m  1. The mean reconstruction error of the temperature field and the reconstruction errors of radiative properties are shown in Fig. 16. As shown in Fig. 15, the reconstruction accuracy of the temperature field and radiative properties decreases with the increase in scattering albedo, especially for the scattering albedo is greater than 0.75. As illustrated in Fig. 16, the reconstruction accuracy of the temperature field and radiative properties decreases with the increase in extinction coefficient. It indicates that the larger the extinction coefficient becomes, the less internal information can be detected, thereby reducing the reconstruction accuracy of the temperature field and radiative properties.

5. Experimental research A non-contact measurement experimental platform of flame temperature field based on the light field imaging has been set up to validate the reconstruction model. The light-field image (L-F image) of ethylene diffusion flame filled with Al2O3 particle has obtained by this experimental platform, and the temperature field and radiative properties of the flame has been reconstructed by using the proposed algorithm. The detailed introduction and parameters of the experimental facility are shown in Ref. [19]. The experimental operating conditions are set as follows. The average particle size of Al2O3 is 30 nm, the flow rate of Al2O3 is 0.1 kg/h, the volumetric flow rate of ethylene is 20 mL/s, the volumetric flow rate of air is 1 L/s, the height of the flame is 150 mm, the radius is 22 mm,

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

55

Fig. 17. The model of simulate light field imaging

of light field imaging, the different directions of radiative intensity are separated by the micro lens and recorded on the pixel units which are covered by the micro lens, so one of the pixel units corresponding to a detection line in the experimental measurements. The direction of detection lines can be determined according to the geometry of the imaging equipment. The acceptance angle of the detector is very small ( o0.015°) [19]. Moreover, the temperature distribution grid is discretized as Nϕ  Nr  Nz ¼1  6  6. The original L-F image obtained by light-field camera and the refocusing image are shown as Fig. 18. The reconstructed 3D temperature distribution in cross sections and vertical section are shown as Fig. 19. From these figures, it can be seen that the temperature of the flame is within the range of 1020 K to 1420 K and similar ranges were also found by others researchers [38,39]. The reconstructed absorption coefficient is 86.41 m  1, and the reconstructed scattering coefficient is 0.02 m  1. Meanwhile, the variation trend of the reconstructed 3D temperature distribution are reasonable, which agrees with the distribution of the light and dark areas of flame (see Figs. 18 and 19). Hence the simultaneous reconstruction method proposed in this manuscript based on the light field imaging technology is practicable for reconstructing the 3D temperature field and radiative properties of flame.

6. Conclusions

Original L-F image

Refocusing image

Fig. 18. The original L-F image (a) and refocusing image (b) of the flame

and the number of the detecting lines is 16,383. In the experiment, the model of simulate light field imaging is established, as shown in Fig. 17. According to the principle

In this research, the multi-dimensional temperature distributions and radiative properties (absorption and scattering coefficients) in absorbing, emitting, and scattering media are reconstructed using GSMFM combined with the hybrid LSQR– SPSO algorithm. A non-contact measurement experimental platform of flame temperature field based on the light field imaging has been set up to validate the reconstruction model and the temperature field and radiative properties of the flame filled with Al2O3 particles has been reconstructed

56

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

Ta/K 1420

150

145 1380 125

120

1340

105

1300

1220 65

Z / mm

Z / mm

90 85

1380 1340

1300 1260

T/K 1420

1260 1220

60

1180

1180

1140

45

1140

30 1100

25

5 -11 -11 R/ 0 0 mm 11 11 mm R/

1060

1060 1020

Cross sections

1100

0 -11 0 11 R / mm

1020

Vertical section

Fig. 19. The reconstructed 3D temperature distribution in cross sections (a) and vertical section (b)

reasonably by using the proposed algorithm. The results indicate that the proposed methodology is a promising approach for recovering realistic temperature distributions. The concluding remarks are summarized as follows. (1) GSMFM is developed as a direct problem-solving method for rapidly simulating the outgoing radiative intensities at the boundary surfaces of participating media. It is proved to have high computing efficiency and accuracy compared with the BMC method. (2) A hybrid LSQR–SPSO algorithm based on the forward GSMFM solution is developed as the inverse problemsolving method for simultaneously reconstructing the multi-dimensional temperature distribution and radiative properties. The axisymmetric temperature distribution and non-axisymmetric bimodal temperature distribution are examined. All of the retrieval results indicate that the temperature distribution and the absorption and scattering coefficients can be retrieved efficiently and accurately using the LSQR– SPSO algorithm combined with GSMFM even with measurement errors. Moreover, the temperature field could be retrieved more accurately than the radiative properties by using the hybrid LSQR–SPSO algorithm. (3) The effects of extinction coefficient and scattering albedo on the accuracy of the estimation are also investigated. The results indicate that the

reconstruction accuracy decreases with the increase of extinction coefficient and the scattering albedo. Further study will focus on simultaneously reconstructing the 3D temperature field and the non-uniform radiative properties field of the absorbing, emitting, and scattering media using the hybrid LSQR–SPSO algorithm. The regularization method will be employed to solve the ill-posed problem for reconstructing the non-uniform temperature distribution and radiative properties.

Acknowledgments The supports of this work by the National Natural Science Foundation of China (No. 51476043, 51576053), and the Major National Scientific Instruments and Equipment Development Special Foundation of China (No. 51327803) are gratefully acknowledged. A very special acknowledgment is also made to the editors and referees who make important comments to improve this paper.

References [1] Liu LH, Tan HP, Yu QZ. Simultaneous identification of temperature profile and wall emissivities in one-dimensional semitransparent

C.-Y. Niu et al. / Journal of Quantitative Spectroscopy & Radiative Transfer 184 (2016) 44–57

[2]

[3]

[4]

[5]

[6]

[7]

[8]

[9] [10]

[11]

[12]

[13]

[14]

[15]

[16]

[17]

[18]

[19]

medium by inverse radiation analysis. Numer Heat Transf Part A 1999;36:511–25. Wang FQ, Tan JY, Wang ZQ. Heat transfer analysis of porous media receiver with different transport and thermophysical models using mixture as feeding gas. Energy Convers Manag 2014;83:159–66. Wu XC, Wu YC, Yang J, Wang ZH, Zhou BW, Grehan G, Cen KF. Modified convolution method to reconstruct particle hologram with an elliptical Gaussian Beam illumination. Opt Express 2013;21: 12803–14. Cui M, Gao XW, Chen HG. Inverse radiation analysis in an absorbing, emitting and non-gray participating medium. Int J Therm Sci 2011;50:898–905. Liu D, Wang F, Yan JH, Huang QX, Chi Y, Cen KF. Inverse radiation problem of temperature field in three-dimensional rectangular enclosure containing inhomogeneous, anisotropically scattering media. Int J Heat Mass Transf 2008;51:3434–41. Cheng Q, Zhang XY, Wang ZC, Zhou HC, Shao S. Simultaneous measurement of three-dimensional temperature distributions and radiative properties based on radiation image processing technology in a gas-fired pilot tubular furnace. Heat Transf Eng 2014;35:770–9. Liu ZJ, Guo Q, Xu L, Ahmad MA, Liu ST. Double image encryption by using iterative random binary encoding in gyrator domains. Opt Express 2010;11:12033–43. Wang F, Yan JH, Cen KF, Huang QX, Liu D, Chi Y, Ni MJ. Simultaneous measurements of two-dimensional temperature and particle concentration distribution from the image of the pulverized coal flame. Fuel 2010;89:202–11. Daun KJ. Infrared species limited data tomography through Tikhonov reconstruction. JQSRT 2010;111:105–15. Qi H, Ruan LM, Shi M, An W, Tan HP. Application of multi-phase particle swarm optimization technique to inverse radiation problem. JQSRT 2008;109:476–93. Zhou HC, Hou YB, Chen DL, Zheng CG. An inverse radiative transfer problem of simultaneously estimating profiles of temperature and radiative parameters from boundary intensity and temperature measurements. JQSRT 2002;74:605–20. Cui M, Gao XW, Chen HG. A new inverse approach for the equivalent gray radiative property of a non-gray medium using a modified zonal method and the complex-variable-differentiation method. JQSRT 2011;112:1336–42. Liu D, Wang F, Cen KF, Yan JH, Huang QX, Chi Y. Noncontact temperature measurement by means of CCD cameras in a participating medium. Opt Lett 2008;33:422–4. Guo Z, Wan SK, August DA, Ying J, Dunn SM, Semmlow JL. Optical imaging of breast tumor through temporal log-slope difference mappings. Comput Biol Med 2006;36:209–23. Akamatsu M, Guo Z. Transient prediction of radiation response in a 3-D scattering-absorbing medium subjected to a collimated short square pulse train. Numer Heat Transf Part A – Appl 2013;63: 327–46. Wang F, Liu D, Cen KF, Yan JH, Huang QX, Chi Y. Efficient inverse radiation analysis of temperature distribution in participating medium based on backward Monte Carlo method. JQSRT 2008;109: 2171–81. Liu D, Yan JH, Wang F. Experimental reconstructions of flame temperature distributions in laboratory-scale and large-scale pulverized-coal fired furnaces by inverse radiation analysis. Fuel 2012;93: 397–403. Liu D, Yan JH, Cen KF. On the treatment of non-optimal regularization parameter influence on temperature distribution reconstruction accuracy in participating medium. Int J Heat Mass Transf 2012;55: 1553–60. Sun J, Xu CL, Zhang B, Hossain MM, Wang SM, Qi H, Tan HP. Threedimensional temperature field measurement of flame using a single light field camera. Opt Express 2016;24:1118–32.

57

[20] McCormick NJ. Inverse radiative transfer problems: a review. Nucl Sci Eng 1992;112:185–98. [21] Liu LH. Simultaneous identification of temperature profile and absorption coefficient in one-dimensional semitransparent medium by inverse radiation analysis. Int Commun Heat Mass Transf 2000;27:635–43. [22] Lou C, Zhou HC. Decoupled reconstruction method for simultaneous estimation of temperatures and radiative properties in a 1-D, gray, participating medium. Numer Heat Transf Part B – Fundam 2007;51: 275–92. [23] Liu D, Yan JH, Wang F, et al. Inverse radiation analysis of simultaneous estimation of temperature field and radiative properties in a two-dimensional participating medium. Int J Heat Mass Transf 2010;53:4474–81. [24] Liu D, Huang QX, Wang F, Chi Y, Cen KF, Yan JH. Simultaneous measurement of three-dimensional soot temperature and volume fraction fields in axisymmetric or asymmetric small unconfined flames with CCD cameras. J Heat Transf 2010;132:1–7. [25] Liu D, Yan JH, Wang F, Huang QX, Chi Y, Cen KF. Simultaneous experimental reconstruction of three-dimensional flame soot temperature and volume fraction distributions. Acta Phys Sin 2011;6: 1–8. [26] Huang Q, Wang F, Yan JH, Chi Y. Simultaneous estimation of the 3-D soot temperature and volume fraction distributions in asymmetric flames using high-speed stereoscopic images. Appl Opt 2012;51(15): 2968–78. [27] Snelling DR, Thomson KA, Smallwood GJ, Gülder ÖL, Weckman EJ, Fraser RA. Spectrally resolved measurement of flame radiation to determine soot temperature and concentration. AIAA J 2002;40: 1789–95. [28] Twynstra MG, Daun KJ, Waslander SL. Line-of-sight-attenuation chemical species tomography through the level set method. JQSRT 2014;143:25–34. [29] Daun KJ, Grauer SJ, Hadwin PJ. Chemical species tomography of turbulent flows: discrete ill-posed and rank deficient problems and the use of prior information. JQSRT 2016;172:58–74. [30] Qi H, Wang DL, Huang XZ, Ruan LM, Zhao H. Development of a general multi-flux method for simulating the radiative intensity in arbitrary direction. J Eng Thermophys 2009;7:1204–6. [31] Paige CC, Saunders MA. LSQR: sparse linear equations and least squares problems. AMC Trans Math 1982;8:195–209. [32] Kennedy J, Eberhart RC. Particle swarm optimization. New York: IEEE; 1995. [33] Qi H, Wang DL, Wang SG, Ruan LM. Inverse transient radiation analysis in one-dimensional non-homogeneous participating slabs using particle swarm optimization algorithms. JQSRT 2011;112: 2507–19. [34] Liu FB. Particle swarm optimization-based algorithms for solving inverse heat conduction problems of estimating surface heat flux. Int J Heat Mass Transf 2012;55:2062–8. [35] Yuan Y, Yi HL, Shuai Y, Liu B, Tan HP. Inverse problem for aerosol particle size distribution using SPSO associated with multilognormal distribution model. Atmos Environ 2011;45:4892–7. [36] Yuan Y, Yi HL, Shuai Y, Wang FQ, Tan HP. Inverse problem for particle size distributions of atmospheric aerosols using stochastic particle swarm optimization. JQSRT 2010;111:2106–14. [37] Zeng JC, Cui ZH. A guaranteed global convergence particle swarm optimizer. J Comput Res Dev 2004;41:1333–8. [38] Lou C, Sun Y, Zhou H. Measurement of temperature and soot concentration in a diffusion flame by image processing. J Eng Thermophys 2010;31(9):1595–8. [39] Ayrancı I, Rodolphe V, Nevin S, Frédéric A, Dany E. Determination of soot temperature, volume fraction and refractive index from flame emission spectrometry. JQSRT 2007;104(2):266–76.