1. 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–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–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.

2. Mathematical formulation and algorithm description

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.

 

Fig. 1 Mathematical formulation of a nonlinear tomography problem.

Download Full Size | PDF

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:

To better elucidate the problem and describe the NIRT algorithm, Eq. (1) is expanded into its element form as shown in Eq. (2):

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:

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:

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.

3. 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(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³, filled with a fluorescence dye solution (a well stirred mixture of ethanol and Rhodamine 6G) to experimentally create a uniform phantom.

 

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

Download Full Size | PDF

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₀ in Eq. (3) was needed to perform the tomographic reconstruction. To obtain I₀, 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₀).

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⁻¹ for this sample.

 

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.

Download Full Size | PDF

Based on the projections obtained on all 7 cameras (each with 800 × 800 pixels, resulting M ≈4.5 × 10⁷ 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³ 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⁶ 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.

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):

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 shot-to-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.

 

Fig. 4 Numerical validation using a uniform phantom. Panel (a): reconstructed distribution on the central plane (i.e., Y = 0) using simulated projections with 4% artificial noises added. Panel (b): reconstructed distribution along 3 lines.

Download Full Size | PDF

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.

 

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

Download Full Size | PDF

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.

4. 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₂) 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₂) jet flow seeded with I₂ vapor. The exit diameter of the jet was 6.35 mm. The I₂ 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₂ 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₂ 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.

 

Fig. 6 A set of example projections on the turbulent flow.

Download Full Size | PDF

Based on the measured projections, 3D distribution of the I₂ vapor concentration in the jet flow was reconstructed using NIRT. These reconstructions were performed in a computational domain of 40 × 40 × 40 mm³ encompassing the entire target flow. The computational domain was discretized into 120 × 120 × 120 voxels (~1.8 × 10⁶), 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₂ concentration, illustrating that the reconstruction captured the overall V-shape of the flow. To further examine the reconstruction, Fig. 7(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₂ 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 Application of NIRT on 3D measurement of I₂ 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.

Download Full Size | PDF

To further elucidate the usefulness of the 3D reconstruction, Figs. 7(c) and 7(d) show a comparison between the 3D reconstructed I₂ 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 shown in Fig. 7(c) was spatially resolved in 3D while that in Fig. 7(d) was line-of-sight-integrated. Such difference was further displayed by comparing the contours of the flow. The white line shown in Fig. 7(c) shows the contour of the flow extracted from the 3D reconstruction, and this contour was overlaid on the line-of-sight-integrated measurements shown in Fig. 7(d).

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 shot-to-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.

 

Fig. 8 Reconstruction error using simulated projections based on turbulent phantoms.

Download Full Size | PDF

 

Fig. 9 The quantitative comparison from several slices and lines between phantom and reconstruction on turbulent V-flow phantom

Download Full Size | PDF

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₂ vapor (i.e., the true answer) across a plane, the second column shows the I₂ 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 shown in Fig. 9(a), the e_R along the second line (5.6%) was significantly smaller than the overall 3D e_R, and the e_R along the third line (10.8%) was significantly larger than the overall 3D e_R. As seen, at a reconstruction accuracy within 8%, the reconstruction was able to capture the overall features of the flow and also most of the detailed features.

Lastly, we conclude this discussion with a discussion of possible approaches that can further improve the reconstruction accuracy of nonlinear tomography problems. As mentioned above, the reconstruction algorithm need to be able to address at least three sources of noises encountered in practice: uncertainty with the camera, the laser system, and the view registration process. The first two sources need to be addressed primarily by experimental equipment and are not discussed further here. The last source can be addressed both by experimental approaches (e.g., design of a new view registration procedure to reduce Δθ) and also by the reconstruction algorithm. One possible approach involves modifying the NIRT algorithm to include the orientations of the cameras (i.e., θ’s) as variables too, so that the view registration errors are accounted for during the reconstruction process. Our preliminary study shows that this approach along could reduce Δθ from = 0.6° to an equivalent of 0.3°, leading to an improved e_R of ~6.6% from ~8% as shown in Fig. 8.

5. Summary

In summary, this work reports the development and validation of a reconstruction algorithm for nonlinear tomography problems, motivated by a range of practical applications involving absorption, scattering, or radiation trapping. Algorithms developed in the past are not readily extendable to such nonlinear problems due to several challenges, such as computation cost, limited view angles available, and the measurement uncertainty. To address these challenges, this work therefore developed a nonlinear iterative reconstruction technique (NIRT). The NIRT algorithm was demonstrated and validated using both experimental tests and numerical simulations. Results show that the NIRT algorithm is able to solve large scale nonlinear tomography problems in 3D effectively with high fidelity.