Development and validation of a reconstruction algorithm for three-dimensional nonlinear tomography problems

This work reports the development and experimental validation of a reconstruction algorithm for three-dimensional (3D) nonlinear tomography problems. Many optical tomography problems encountered in practice are nonlinear, for example, due to significant absorption, multiple-scattering, or radiation trapping. Past research efforts have predominately focused on reconstruction algorithms for linear problems, and these algorithms are not readily extendable to nonlinear problems due to several challenges. These challenges include the computational cost caused by the nonlinearity (which was compounded by the large scale of the problems when they are 3D), the limited view angles available in many practical applications, and the measurement uncertainty. A new algorithm was therefore developed to overcome these challenges. The algorithm was validated both numerically and experimentally, and was demonstrated to be able to solve a range of nonlinear tomography problems with significantly enhanced efficiency and accuracy compared to existing algorithms. ©2016 Optical Society of America OCIS codes: (110.6960) Tomography; (100.6890) Three-dimensional image processing. References and links 1. S. R. Arridge, “Optical tomography in medical imaging,” Inverse Probl. 15(2), R41–R93 (1999). 2. W. Yang and L. Peng, “Image reconstruction algorithms for electrical capacitance tomography,” Meas. Sci. Technol. 14(1), R1–R13 (2003). 3. R. Santoro, H. Semerjian, P. Emmerman, and R. Goulard, “Optical tomography for flow field diagnostics,” Int. J. Heat Mass Transfer 24(7), 1139–1150 (1981). 4. Y. Ishino and N. Ohiwa, “Three-dimensional computerized tomographic reconstruction of instantaneous distribution of chemiluminescence of a turbulent premixed flame,” JSME Int. J. Ser. B 48(1), 34–40 (2005). 5. J. Floyd, P. Geipel, and A. M. Kempf, “Computed tomography of chemiluminescence (CTC): Instantaneous 3D measurements and Phantom studies of a turbulent opposed jet flame,” Combust. Flame 158(2), 376–391 (2011). 6. J. Floyd and A. Kempf, “Computed tomography of chemiluminescence (CTC): high resolution and instantaneous 3-D measurements of a matrix burner,” Proc. Combust. Inst. 33(1), 751–758 (2011). 7. L. Ma, Q. Lei, Y. Wu, W. Xu, T. M. Ombrello, and C. D. Carter, “From ignition to stable combustion in a cavity flameholder studied via 3D tomographic chemiluminescence at 20 kHz,” Combust. Flame 165, 1–10 (2016). 8. W. Cai, X. Li, F. Li, and L. Ma, “Numerical and experimental validation of a three-dimensional combustion diagnostic based on tomographic chemiluminescence,” Opt. Express 21(6), 7050–7064 (2013). 9. Q. Lei, Y. Wu, H. Xiao, and L. Ma, “Analysis of four-dimensional Mie imaging using fiber-based endoscopes,” Appl. Opt. 53(28), 6389–6398 (2014). 10. Y. Wu, Q. Lei, and L. Ma, “Experimental demonstration of 4D imaging in two-phase flows based on computed tomography at 5 kHz,” Appl. Opt. 53(24), 5547–5553 (2014). 11. F. Scarano, “Tomographic PIV: principles and practice,” Meas. Sci. Technol. 24(1), 012001 (2013). 12. C. Atkinson and J. Soria, “An efficient simultaneous reconstruction technique for tomographic particle image velocimetry,” Exp. Fluids 47(4-5), 553–568 (2009). 13. M. M. Hossain, G. Lu, and Y. Yan, “Optical fiber imaging based tomographic reconstruction of burner flames,” IEEE Trans. Instrum. Meas. 61(5), 1417–1425 (2012). 14. W. Cai, D. J. Ewing, and L. Ma, “Investigation of temperature parallel simulated annealing for optimizing continuous functions with application to hyperspectral tomography,” Appl. Math. Comput. 217(12), 5754–5767 (2011). 15. K. Belkebir, P. C. Chaumet, and A. Sentenac, “Influence of multiple scattering on three-dimensional imaging with optical diffraction tomography,” J. Opt. Soc. Am. A 23(3), 586–595 (2006). #266609 http://dx.doi.org/10.1364/OE.24.015912 Journal © 2016 Received 19 May 2016; revised 30 Jun 2016; accepted 30 Jun 2016; published 6 Jul 2016 Vol. 24, No. 14 | 11 Jul 2016 | OPTICS EXPRESS 15912


Introduction
An incredibly wide range of practical problems can be mathematically formulated as a tomography problem.The goal of tomography is to reconstruct a spatially-resolved distribution of a certain property of a target based on a set of line-of-sight integrated measurements of the target measured from different perspectives (and these measurements are called projections).Different applications involve different target and physical processes, and examples range from the medical imaging of human body using optical transmission measurements [1], examination of the internal structure of solids or liquids using electrical capacitance measurements [2], and optical imaging of a chemical species in fluid flows using absorption or emission measurements [3,4].Regardless of such wide range of applications and the corresponding physical processes involved, the mathematical background involved remains essentially the same.Hence, the rest of paper will be developed under the context of optical imaging of fluid flows, even though the mathematics and algorithm can be extended to other applications straightforwardly.
Mathematically, tomography problems can be broadly divided into linear and nonlinear problems.For a linear tomography problem, the projections depend linearly on the property to be reconstructed.In the area of optical imaging of fluid flows, an example of a linear problem involves the imaging of the concentration of a radical (such as CH*) using tomographic chemiluminescence when radiation trapping is negligible [4][5][6][7][8].Because of the negligible radiation trapping, the projections (i.e., line-of-sight-integrated signals emitted by CH*) are proportional to the sought property (i.e., the concentration of CH*).When radiation trapping is not negligible, the projections will depend nonlinearly on the sought property and the problem becomes nonlinear.Example of other linear problems include tomographic Mie scattering when multiple-scattering is negligible [9,10] and tomographic particle image velocimetry [11,12] (which is essentially tomographic Mie scattering with negligible multiple-scattering). Extensive research efforts have invested in the development of reconstruction algorithms for solving linear tomography problems.A few widely established algorithms include the algebraic reconstruction technique (ART) [5][6][7], multiplicative algebraic reconstruction technique (MART) [11], filtered backprojection (FBP) [13] and reconstruction technique based on minimization [8,14].For a nonlinear tomography problem, the projections depend nonlinearly on the property to be reconstructed.Besides the example of tomographic chemiluminescence when radiation trapping is not negligible mentioned above, other nonlinear problems include tomographic scattering measurements with appreciably multiple-scattering [15,16], and tomographic laser induced fluorescence (LIF) measurements with significant absorption [17] or amplified spontaneous emission [18].
Nonlinear problems pose several unique challenges, and reconstruction algorithms developed for linear problems are not easily extendable to solve nonlinear problems.Algorithms proposed in the past for nonlinear tomography problems could be broadly divided into two categories.Algorithms in the first category are based on the central slice theorem [19], an specific algorithms include the convolution methods [20,21] and filtered backprojection [22,23].Algorithms in this category require a large number of projections from a wide range of perspectives (hundreds and more).However, in many practical problems (such as flow imaging), only a limited number of projections is available, typically ranging from 2 [24] to about 50 [4].It is practically difficult or infeasible to obtain a large number of projections from a wide angular range due to the transiency of the target, the optical access available, and the cost of equipment.Algorithms in the second category are based on the conversion of the tomography problem into a minimization problem, for example, to minimize a cost function defined as the difference between the measured and calculated projections.Algorithms in this category enjoy the flexibility of being applicable to both linear and nonlinear problems.However, in practice, the problem may involve a large number of variables (e.g., on the order of millions for the problems examined in this work) and the projections are measured with certain uncertainty.As a result, it is not trivial to find a minimization technique that can minimize the cost function with efficiency, robustness, and accuracy, either for linear or nonlinear problems [25].Minimization techniques investigated in the past included approaches based on the gradients of the cost function, such as the conjugate gradient technique [15,16,26,27], Gauss-Newton technique [28], and Newton-Kantorovich technique [29,30].Unfortunately, practical problems often feature many confusion local minima (even more so for nonlinear problems than linear problems) and gradient-based techniques can only converge to one of these many local minima within the measurement uncertainty [31].To overcome this limitation, global optimization techniques such as simulating annealing (SA) [25,32,33] and genetic algorithm (GA) [34] have been adopted to solve the tomography problems.These algorithms, at least in principle, can find the global minimum or a solution sufficiently close to it [35].However, the computational cost of those techniques is dramatically higher than those based on gradients and can become prohibitive for problems with practical scale.
Based on the above understanding of past work, this work describes a new approach, termed NIRT (Nonlinear Iterative Reconstruction Technique), to solve nonlinear tomography problems.The NIRT approach solves the nonlinear problem iteratively in a manner analogous to the solution of linear tomography problem using ART.The NIRT algorithm is demonstrated to solve nonlinear tomography problems while retaining the advantages of ART for linear problems.For example, the NIRT algorithm was able to solve nonlinear problem with limited number of projections, with robustness in the presence of uncertain projections, and with good computational efficiency.In the rest of this paper, we first detail the mathematical formulation of the NIRT algorithm, and then report its numerical and experimental validation.This section introduces the mathematical formulation of nonlinear tomography and describes the NIRT algorithm.Figure 1 schematically illustrates the mathematical formulation.As can be seen, the 3D distribution of the sought property is denoted as C (e.g., the concentration distribution of a chemical species or particulates) and discretized into voxels under a Cartesian coordinate system.An imaging system collects the signal emitted by the target property to form an image on the camera chip (defined as a projection and denoted as P).As an example, the dashed red and brown lines illustrate that the signals emitted by two different voxels arrive on two sets of pixels on the camera chip and that these two sets of pixels can overlap.The projection P depends on two factors: the parameters of the imaging system and the modulation (attenuation or magnification) of the signal by the target property itself.The parameters of the imaging system include the distance and orientation of the imaging system, specified by r (distance), θ (azimuth angle), and φ (inclination angle).The effects of all these parameters are reflected in the point spread function (PSF) of the imaging system [36,37].The modulation of the target signal could be caused by a variety of physical processes depending on the specific signal generation and propagation mechanisms.As a simple example shown in Fig. 1, when the signal generation involves a laser (e.g., in LIF or Mie scattering measurements, respectively), the signal could be attenuated by the absorption or scattering of the illumination laser or magnified by amplified spontaneous emission.When the signal generation does not involve a laser (e.g., in chemiluminescence measurements), the signal could be attenuated by radiation-trapping.

Mathematical formulation and algorithm description
Regardless of the specific signal generation and modulation mechanism, the relationship between the measured projection P and the sought property can be mathematically capsulated into the following equation: where P presents the projection vector formed by organizing the measured projections pixel by pixel, PSF the point spread function in matrix form, C the vector formed by organizing the discretized sought property voxel by voxel, and lastly Г a modulation function representing the signal modulation mechanism.In Eq. ( 1), the operator '•' represents the Hadamard product and '•' matrix and vector product.Note that the Hadamard product is associative but the matrix product is not.This work uses bold letters to symbolize vectors or matrices formed by discretizing their corresponding continuous functions.The tomographic problem is to solve for C with P measured at different locations and orientations.Based on Eq. ( 1), it can be seen that when Г does not depend on C or depends linearly on C, the tomographic problem is a linear problem.This work focuses on nonlinear problems where Г depends nonlinearly on C.
To better elucidate the problem and describe the NIRT algorithm, Eq. ( 1) is expanded into its element form as shown in Eq. ( 2): where m and M are the index and total number of pixels on the projection, respectively; P m is measured projection on the m th pixel; n and N are the index and total number of voxels in the measurement domain, respectively; PSF m,n is the value of PSF matrix for the n th voxel on the m th pixel, and C n and Г n are the values of the sought property and the nonlinear modulation function for the n th voxel, respectively.Since the PSF matrix does not depend on C, established ART algorithm [38] can be applied to solve for Г•C voxel by voxel.However, in general, Г•C cannot be decoupled by ART to obtain C when Г depends on C nonlinearly.As an example, Г assumes the following exponential form according to the Beer-Lambert relationship for a variety of signal generation processes (such as absorption, scattering, or radiation trapping): where α represents the attenuation coefficient (either absorption or scattering coefficient), l the pathlength propagated by the illumination laser and/or the signal photons before reaching the imaging system, C l the value of C along l, and I 0,n the incident laser intensity on the n th voxel.Here, Γ depends not only on C nonlinearly but also on l, both may vary from voxel to voxel.The established ART algorithm is unable to address such dependence in general.
The NIRT algorithm was developed to address these challenges and solve the nonlinear tomography problem described above.The first step of the development was to cast Eq. ( 1) into a different form as shown below: by converting Γ from a vector to a diagonal matrix (i.e., diag(Γ)) with the elements of vector Г on the main diagonal.The second step was to define a new point spread function, NPSF, by combining the first two terms of right hand in Eq. ( 4).The specific element form of NPSF is shown below: , Conceptually, NPSF represents the nonlinear point spread function that takes both the geometric PSF and the nonlinear modulation Г into account.Now Eq. ( 1) can be rewritten into: The third step was to solve Eq. ( 6) iteratively according to the following scheme: where q represents the index of the iteration step, q m D the difference between the measured (P m ) and calculated projection ( q m P ) at the m th pixel in the q th iteration, β a relaxation factor set to be 0.01 in this work, and || || q NPSF the L2 norm of NPSF in the q th iteration.As can be seen from the Eqs.( 7) and 8, the NIRT algorithm accounts for the modulation of the signal in each step of the iteration to address the nonlinearity (at the cost of increased computational cost compared to the linear ART algorithm, since the NPSF matrix needs to be updated in every iteration).To address the signal modulation's dependence on l (i.e., the propagation pathlength), a ray-tracing module based on the Monte Carlo method [18] was employed to track the signal modulation on a voxel by voxel bases.This ray-tracing module is an important component of the NIRT algorithm, and it enables the NIRT algorithm to be applicable to a wider range of tomographic problems compared to past algorithms.For instance, past algorithms typically were developed under the so-called parallel beam or fan beam assumption.These assumptions have limited range of validity in practice, primarily in two aspects.First, these assumptions become invalid in practice when the signals emitted from different voxels overlap (as illustrated in Fig. 1) due to focusing and zooming.Second, algorithms based on these assumptions are typically limited to two-dimensional (2D) problems or 3D problems that can be decouple into a series of 2D problems, while many practical problems are truly 3D due to signal overlapping and pathlength dependent signal modulation, and cannot be decoupled into a series of 2D problems.
Lastly, the iteration terminates when the change in C during two consecutive iterations is less than a preset value.More specifically, the following termination criterion as described in [39] was employed in this work: where Δ is a small number that controls the speed of convergence.The value of Δ was empirically suggested to be in a range of 0.0001% to 0.1% from previous work [36,39], and was preset to be 0.1% in this work.The results from this work showed that the reconstruction accuracy was insensitive to the choice of Δ (the reconstruction accuracy changes within 0.5% when Δ varied between 0.001% to 0.1%).Note that there is no extra constraint placed on C by the algorithm besides that it should be positive.In summary, the development of the NIRT algorithm relied on two key enabling techniques to overcome the challenges of the nonlinear tomography problems.The first enabling technique involved the derivation and definition of NPSF, the nonlinear point spread function, which allowed the development of the iterative scheme to decouple the nonlinear problem.The second enabling technique involved the application of a ray-tracing module to account for dependence of signal modulation on pathlength.The rest of this paper will describe the experimental and numerical validation of the algorithm, and then the demonstration of the algorithm in an application involving volumetric LIF measurement in turbulent flows.

Validation of algorithm
The NIRT algorithm has been first extensively validated numerically using a range of phantoms.Here, we will focus on the description of the numerical validation on a phantom for which an accompanying experimental validation was performed, and only briefly discuss the numerical validation on other phantoms.This validation involved a phantom with a uniform distribution of the target species, which can be experimentally realized with excellent control and accuracy as shown in Fig. 2(a).Figure 2  When illuminated volumetrically by a laser at a wavelength of 532 nm (generated by a pulsed Nd:YAG and expanded by a series of lenses as shown), the dye solution absorbs the illumination laser as it propagates through the solution and emits LIF photons, resulting in a nonlinear problem as described above.The dye solution emits LIF signals volumetrically, which were captured by a total of 7 cameras (6 Photron SA4s and 1 SA6) from different orientations.The goal of this validation was to compare the reconstructed distribution of the dye concentration to the uniform distribution known a priori.
To facilitate the following discussion, a right-handed Cartesian coordinate system was defined as shown in the figure: the origin was defined as the center point of the cubical cell, the X axis was defined along the propagation direction of the illumination laser (which propagated perpendicularly into the dye cell), and the Z axis was defined to be out of the paper as shown in Fig. 2(a).In the experiment, all 7 cameras were aligned in the X-Y plane and therefore their orientations were completely specified by θ, defined as the angle formed by the optical axis of a given camera relative to the positive X direction as shown.The focal length and f-number of the lenses used on all cameras were 105 mm and 2.8, respectively.Each lens was equipped with a 532 ± 10 nm OD4 notch filter to block any interference due to the scattering or reflection of the illumination laser.Prior to any measurement, a calibration target was used in view registration program to determine the orientations of cameras [40].The orientations of camera 1 through 7 were determined to be θ = 90.0°,127.5°, 213.8°, 240.2°, 272.2°, 316.9° and 48.9°, respectively, with an accuracy estimated to be within 0.6°.
The profile of the illumination laser pulses (i.e., I 0 in Eq. ( 3) was needed to perform the tomographic reconstruction.To obtain I 0 , another dye cell, with a length and height both of 50 mm and a thickness of 0.5 mm, was fabricated.This cell was used to hold a thin layer of the same dye solution, and was then placed at X = 25 mm, i.e., where the laser entered the target measurement volume.The thin cell was illuminated by the laser, and an image of the LIF signal emitted was captured as shown in Fig. 2(b).Due to the thinness of the dye solution in the cell, the integration effects in the dye were neglected and the image shown in Fig. 2(b) was used to represent the intensity profile of the illumination laser (i.e., I 0 ).To confirm the absorption, a camera was placed to capture the projection of the LIF signal emitted by the dye cell from the top view, as shown in Fig. 3(a).If the absorption was negligible, the LIF signal captured from this view should be constant along the X direction (i.e., the propagation direction of the laser) due to the uniformity of the dye concentration.Figure 3(b) shows the intensity of the signal along the X direction at three Y locations (Y = −5, 0 and 5 mm, illustrated by the dashed lines in panel a).As seen in Fig. 3(b), signal decreased appreciably as the laser propagated into the dye solution, confirming that the absorption was not negligible.A curve fitting on the data shown in Fig. 3(b) also confirmed that the signal decreased exponentially according to Beer-Lambert relationship as shown in Eq. ( 3).The absorption coefficient, α, was estimated to be around 0.006 mm −1 for this sample.
Based on the projections obtained on all 7 cameras (each with 800 × 800 pixels, resulting M ≈4.5 × 10 7 as shown in Eq. ( 2), a tomography reconstruction was performed using the NIRT algorithm to obtain the 3D distribution of dye concentration (which should be ideally a uniform distribution because the dye solution was well mixed).In this volumetric LIF problem, the two-level LIF model was used and the fluorescence yield was set as a constant value for all particles in the concentration field (i.e., the goal is essentially to measure a relative concentration field, which is what many LIF measurements are aimed at).The computational domain was taken to be a volume of 40 × 40 × 40 mm 3 center around the origin, smaller than the volume of the dye cell to exclude non-ideal effect around the edges and corners of the cell.The computational domain was then discretized into 120 × 120 × 120 voxels (i.e., N ≈1.8 × 10 6 as shown in Eq. ( 2), resulting in a nominal resolution of 0.33 mm in all three directions.The nominal resolution was calculated by dividing the dimension of the measurement volume with the discretization (i.e., 40 mm/120 = 0.33 mm).This nominal resolution represents the highest possible resolution that can be obtained.The actual resolution will be worse than the nominal resolution due to the factors such as measurement and reconstruction uncertainty.Our previous work [17] have quantified the spatial resolution by comparing the volumetric LIF with a planar LIF (PLIF) measurement on a turbulent flow.Based on the comparison, the actual spatial resolution was estimated to be 0.71 mm.Also note that M is larger than N in this work, resulting in an over-determined systems (i.e., more equations than unknowns).In practice, we always prefer to design the experiments to make the problem over-determined.The number of equations is controlled by the number cameras and their field of view, and the number of unknowns by the discretization scheme.3(c).As seen from Fig. 3(d), the reconstructed concentration agreed well with the expected uniform distribution (shown as dotted lines).To quantify the reconstruction accuracy, an average reconstruction error (e R ) was defined according to Eq. ( 10): represent the reconstructed and true concentration distribution on the n th voxel, respectively.For results shown in Fig. 3(d), e R was on the order of 4%.
Numerical validation was performed parallel with the above experiments, and such numerical validation also provided further insights to the e R calculated above.For the numerical validation, a uniform phantom was created to simulate the distribution of the dye in the cell.Projections from the same seven orientations as used in the experiments were calculated using the laser profile shown in Fig. 2(b).To simulate the experimental uncertainty, a total of 4% of noise was artificially added to the projections.This noise level was chosen based on an analysis of the experiments.Three primary sources of uncertainty were identified in the experiments: the uncertainty with the camera (such as background noises, camera nonlinearity and/or nonuniformity et al., estimated to be about 1%), the shotto-shot variation of the laser system (estimated to be about 1.1%), and the accuracy of our view registration process.As aforementioned, the accuracy of the view registration was estimated to be within 0.6°, and such an uncertainty translated into about 1.6% uncertainty in the projection.The overall uncertainty from all three sources was therefore estimated to be around 4%. Gaussian noise was added to the numerical projection data to simulate the background noise and the shot-to-shot variation.As for the view registration noise, the noisy projections were simulated using the angle sets with measurement uncertainty.Based on the simulated projections, reconstructions were performed and the results summarized in Fig. 4. Figure 4(a) shows the reconstructed concentration at a slice corresponding to Y = 0 mm using the simulated projections.Figure 4(b) shows the reconstructions along three lines in comparison to the phantom (shown in black dash lines), with e R calculated to be 4.01%, 4.03% and 4.53%, respectively.These results were in good agreement with those obtained experimentally, supporting the accuracy of the algorithm and the analysis of the experimental uncertainty.Lastly, to further elucidate the difficulties associated with nonlinear tomography problems, numerical studies were performed to compare the accuracy and computational cost of three reconstruction techniques: NIRT, ART, and a minimization technique.When the ART algorithm was used, it was used to solve the problem as if it were a linear problem (i.e., the absorption is neglected).In the minimization technique, the tomography problem was converted into a problem to find a distribution that can reproduce the projections with minimum difference.Due to the complexity of this minimization process, a stochastic minimization algorithm based on simulated annealing (SA) [33] was used.Figure 5 compares the e R and computing time when these techniques were applied on the same uniform phantom.All simulations were performed with 4% artificial noise added to the projections, and all computations were performed on a same computational workstation with 16-Core Intel Xeon 2.6GHz CPU (only 1 core is used in all calculations) and 512GB memory.As can be seen, both NIRT and SA achieved the same level of accuracy, indicating that their ability to converge to the globally optimal solution for complicated tomographic problems is comparable.Both NIRT and the minimization technique SA yielded good reconstruction accuracy with e R below 5%, illustrating the effectiveness of NIRT method and also the flexibility of the minimization technique and the prowess of SA to minimize large scale and complicated functions.In contrast, ART failed to solve the problem correctly and resulted in an e R ≈60% due to its inability to address the nonlinearity in the problem.
In terms of the computational cost, the computing time of NIRT (~7,200 seconds) was on the same order of magnitude as that of ART (~5,760 second) because the NPSF needed to be updated in each iteration in NIRT and required extra computation time.The SA algorithm, in contrast, cost approximately 30 × more computing time than NIRT (~219,600 seconds).In conclusion, the NIRT algorithm was capable of performing 3D nonlinear tomography with reconstruction fidelity as methods based on minimization techniques, while with computational cost on the same order as ART.

Demonstration on a turbulent flow
After the above experimental and numerical validations, this section reports an application of the NIRT algorithm on the 3D measurement of iodine (I 2 ) vapor using a technique we codenamed VLIF (volumetric LIF).The experimental setup was similar to that shown in Fig. 2(a), with two major differences.The first major difference was that the controlled cell was replaced by a turbulent nitrogen (N 2 ) jet flow seeded with I 2 vapor.The exit diameter of the jet was 6.35 mm.The I 2 vapor was introduced into the flow by heating solid iodine crystals in a water bath to a temperature of 200 °F.The mole concentration of the I 2 vapor in the target flow was estimated to be 4%.A rod (with a diameter of 3.18 mm) was placed at the exit of the jet, to increase the complexity of the flow structures (and also to create an easily recognizable pattern of the flow to facilitate the discussion of the results).The second difference was that the 10 Hz Nd:YAG laser was replaced by a 5 kHz pulsed laser at a wavelength of 527 nm (Photonics Industries DM20-527).The pulse duration of the 5 kHz laser was ~300 ns and the pulse energy was 24 mJ.Similar to the dye cell measurement, the laser pulses generated were expanded to illuminate the flow volumetrically, exciting the seeded I 2 molecules to generate LIF photons volumetrically.Base on this new setup, Fig. 6 shows a set of example projections of the target flow captured by camera 1 through 7, and also an example 2D laser profile measured using the thin dye cell mentioned before.All images shown here were single-shot measurements.The projections displayed a turbulent jet flow with an overall V shape generated by the rod placed at the exit of the jet.All projections had a resolution of 800 × 800 pixels, and each pixel corresponded to a physical dimension of 0.05 × 0.05 mm.The overall size jet flow was about 40 mm in each direction.The absorption of the illumination laser by the seeded iodine was estimated to be ~5% as it propagated through the measurement volume, not as large as the absorption observed in the dye cells but large enough to preclude the use of linear algorithms such as ART.
Based on the measured projections, 3D distribution of the I 2 vapor concentration in the jet flow was reconstructed using NIRT.These reconstructions were performed in a computational domain of 40 × 40 × 40 mm 3 encompassing the entire target flow.The computational domain was discretized into 120 × 120 × 120 voxels (~1.8 × 10 6 ), resulting in a nominal spatial resolution of 0.33 mm. Figure 7 shows an example reconstruction obtained using the projections shown in Fig. 6. Figure 7(a) shows a 3D rendering of the relative I 2 concentration, illustrating that the reconstruction captured the overall V-shape of the flow.To further examine the reconstruction, Fig. 7  Similar to the validation studies described in Section 3, parallel numerical simulations were performed for these measurements conducted in turbulent flows to provide further insights.In these simulations, a phantom distribution resembling that target turbulent flow was created.Based on the phantom, projections at the same orientations used in the experiments were computed with artificial noises added to simulate measurement uncertainty.Then these simulated projections were used in the NIRT algorithm.As aforementioned, the noises considered included the uncertainty with camera (estimated to be about 1%), the shotto-shot variation of the laser system (estimated to be about 1.1%), and the accuracy of our view registration process.Figure 8 examines the last noise source more closely by showing the e R obtained under various level of uncertainty caused by view registration error.As seen, at our current level of view registration accuracy of Δθ = 0.6°, the error in the projection (e P, VR ) caused by Δθ was about 2.7% (which was larger than that for the simpler uniform distribution examined in Section 3).In general, our results suggest that the impact of Δθ on e P, VR increases as the complexity of the distribution increases, motivating the improvement of the view registration process as to be further elaborated at the end of this section.With an e P, VR of 2.7% and a combined uncertainty of 2.1% due to the background noise registered on the camera and the shot-to-shot variation, a total of 4.8% artificial noise was added to the computed projections to simulate the experiments.The e R obtained was ~8% as shown in Fig. 8, significantly larger than that obtained for the simpler uniform phantom discussed in Section 3 due to the increased complexity in the target distribution.Before further discussion on such reconstruction accuracy and on possible ways to improve it, we wanted to first point out that a reconstruction with an e R of 8% is actually quite accurate already and resolve most of the detailed features as illustrated in Fig. 9. Figure 9 shows a more detailed comparison between the turbulent phantom and the reconstruction using simulated projections (with 4.8% artificial noise added) across several planes and along several lines.In each panel, the first column shows phantom distribution of I 2 vapor (i.e., the true answer) across a plane, the second column shows the I 2 distribution reconstructed by NIRT across the same plane, and the third column compares the phantom and the reconstruction along three lines (the dashed line).These three lines were picked so that the e R along the first line (8.9%) was about the same as the overall 3D e R of 8% (as simulations.Results show that the NIRT algorithm is able to solve large scale nonlinear tomography problems in 3D effectively with high fidelity.
(a)  shows the setup of the validation experiment from the top view.The setup used a cubical cell, with a dimension of 50 × 50 × 50 mm 3 , filled with a fluorescence dye solution (a well stirred mixture of ethanol and Rhodamine 6G) to experimentally create a uniform phantom.

Fig. 2 .
Fig. 2. Panel (a): experimental setup from the top view.Panel (b): the spatial profile of illumination laser intensity at X = 25 mm.

Fig. 3 .
Fig. 3. Experimental results obtained from the controlled dye cell.Panel (a): projection of the volumetric LIF signal from the top view.Panel (b): intensity of signal along the laser propagation direction at three different Y locations (Y = −5, 0 and 5 mm).Panel (c): the reconstructed concentration distribution at three difference planes (Z = 10, 20 and 30 mm).Panel (d): reconstructed dye concentration along three lines.

Figure 3 (
Figure 3(c) shows the result of the reconstruction on three planes (corresponding to Z = 10, 20 and 30 mm).As seen, the reconstruction shows the expected uniform distribution.To examine the reconstruction closely and quantitatively, Fig. 3(d) shows the reconstructed concentration along three lines shown by the dashed line in Fig.3(c).As seen from Fig.3(d), the reconstructed concentration agreed well with the expected uniform distribution (shown as dotted lines).To quantify the reconstruction accuracy, an average reconstruction error (e R ) was defined according to Eq. (10):

Fig. 5 .
Fig. 5. Comparison of reconstruction accuracy and computational time among different algorithms.

Fig. 6 .
Fig. 6.A set of example projections on the turbulent flow.
(b) shows three planar slices of the reconstruction at three different Y locations (Y = −5, 0 and 5 mm).As seen from Fig. 7(b), the 3D reconstruction also captured the flow feature at a more detailed level.For example, the slice at Y = 0 corresponded to the central plane of the jet flow, and therefore the overall I 2 concentration on this slice should be higher than the other two slices due to less entrainment, as clearly captured by the 3D reconstruction shown in Fig. 7(b).

Fig. 7 .
Fig. 7. Application of NIRT on 3D measurement of I 2 concentration in turbulent flows using VLIF.Panel (a): 3D rendering of the reconstruction.Panel (b): three planar slices of the 3D reconstruction.Panel(c): a planar slice taken out of 3D reconstruction at Y = 5 mm.Panel (d): projection measured by camera 5. To further elucidate the usefulness of the 3D reconstruction, Figs.7(c) and 7(d) show a comparison between the 3D reconstructed I 2 distribution against a measured projection.Figure 7(c) shows a planar slice of the 3D reconstruction at Y = 5 mm, and Fig. 7(d) shows the projection measured by camera 5 (which was placed in the same view orientation as Fig. 7(c)).The images shown in Fig. 7(c) and 7(d) differ apparently, due to the fact that the result