The Effect of Acoustic Impedance on Subsurface Absorber Geometry Reconstruction using 1D Frequency-Domain Photoacoustics.

This paper considers the effect of an impedance mismatch between the absorber and its surroundings on the aborber reconstructions from the photoacoustic signal profile, in particular when a non-delta input pulse is used. A transfer function approach is taken, demonstrating in the case of impedance mismatch how the total response can be modeled using the sum of the mismatch-free response and its time-delayed, time-reversed replicas, which may or may not overlap. It is shown how this approach can be exploited to accommodate the effects of non-delta pulses and/or pulse-equivalent waveforms such as linear-frequency-modulated (LFM) chirps, and impedance mismatches in any inversion algorithms, even in the presence of large reflection coefficients. As a consequence, for simple-absorber reconstruction algorithms that assume impulses or 'short enough' pulses, the compressive portion of the measured response may be used in reconstruction formulas that do not model the impedance mismatch, regardless of the size of the mismatch. For longer-duration input waveforms, it is demonstrated how existing reconstruction methods can be successfully adapted to include the effect of the impedance mismatch. Simulations are used to illustrate these ideas. The gained physical insight into how components of the generated pressure wave carry absorber information is then exploited for signal inversion and absorber reconstruction in the frequency domain when multi-frequency modulation chirps are used for photoacoustic radar pressure measurements. The foundational theoretical developments ultimately address impendance mismatch issues germane to the major photoacoustic frequency-domain imaging modality to-date, which is the photoacoustic radar.

Introduction

Photoacoustics is the generation of acoustic waves as a consequence of the absorption of light energy by an absorbing material and the subsequent thermoelastic expansion of the material. This effect has shown great promise for biomedical imaging applications. There has been much interest in PA applications because the method combines the spatial resolution of ultrasound with the contrast of optical absorption for deep imaging in biological tissues [1], [2], [3].

Diebold [4] gave a concise explanation of the governing equation for the pressure that results from launching a photoacoustic wave. Light pulses are delivered into biological tissue, which then absorbs the light and converts it to heat, generating an initial pressure rise due to thermoelastic expansion. The initial pressure gives rise to a photoacoustic signal, which is then detected by a transducer. The governing equation is given bywhere β is the thermal expansion coefficient, c is the speed of sound, C is the specific heat, H is the energy per unit volume and time deposited by the optical radiation beam, and is the pressure of the acoustic wave, a function of space and time. The heating function H is assumed to be the result of a chromophore absorber with an optical absorption coefficient μ (dimensions of inverse length) that is heated by an optical pulse with fluence (energy per unit area) of F. As is common, it is assumed that H is a separable function of space and time so that . Thus, is a function of space that describes the geometry of the absorber. Similarly, I(t) is a function that describes the time dependence of the incident optical pulse.

The impulse response of Eq. (1) has been extensively studied for different geometries in free space with the assumption of a uniform speed of sound, implying the same characteristic impedances for an absorber and its surroundings [5]. In fact, many optoacoustic inversion schemes are derived by assuming delta pulse excitation and also by assuming that the imaged sample is acoustically homogeneous and perfectly matched to the coupling medium. For example, a radon-transform-like approach was used in [6], iterative approaches were used in [7] and [8], a Green's function approach was used in [9], a time-reversal approach was used in [10], and reconstructions from spherical and circular mean data were considered in [11]. Additionally, in recent years, acoustic heterogeneity of the medium surrounding the absorber has been incorporated into modeling approaches. Dean-Ben et al. considered various complex weighting functions to account for acoustic heterogeneities by assuming some a priori information on the distribution of the acoustic heterogeneities in the surrounding medium [12], [13], [14]. Other authors have considered acoustic heterogeneities in time-reversal schemes, even though the heterogeneity may destroy the assumption of a finite time response as waves may be trapped in the medium owing to acoustic reverberations [15], [16].

Fewer efforts have been made to consider the case where the absorber itself possesses different acoustic impedance from its surroundings and its implications on the shape of the detected acoustic pulse. This is particularly of interest when the probing input pulse does not satisfy the delta function distribution commonly assumed for photoacoustic tomography. Diebold et al. discussed the case of homogeneous light absorption with a sinusoidal heating function [17] and pulsed radiation [18]. Karabutov et al. developed a theory for homogeneous, microhomogeneous and microinhomogeneous media for pulsed later photoacoustics [19]. Also, Herrerías-Azcué et al. [20] considered the 1D case of an absorber with Lambert Beer absorption and a delta function impulse heating function. A semi-infinite 1D Beer absorption case was considered in [21] and analytical solutions under the assumption of homogeneous light absorption have been reported for a sphere in [22]. In a similar vein to this paper, Anastasio et al. [23] considered the effect of heterogeneous acoustic properties of the object in thermoacoustic tomography and heuristically introduced the idea of using half-time reconstruction in order to mitigate the effect of the absorber's acoustic heterogeneity. It shall be shown herein that those heuristic ideas have indeed a valid theoretical basis.

In this paper, the effect of an impedance mismatch between the absorber and its surroundings on absorber reconstruction is examined in the case of frequency-domain photoacoustic excitation waveforms which exhibit narrower spectral bandwidths than short (ns) laser pulses. A plane, finitely-thick optically absorbing sample surrounded by a transparent fluid of different acoustic impedance is considered, implying a photoacoustic equation in one dimension. Analytical pressure response expressions are given for the general case, independent of the specific form of the heating function and the depth-dependent absorption coefficient. The implications of the impedance mismatch on the measured reflection and transmission pressure response are examined, with particular emphasis on how the presence of the impedance mismatch, in conjunction with a non-delta function (long) input pulse or pulse-equivalent waveforms, combine to alter the shape of the produced photoacoustic pressure pulse. In particular, we theoretically show for the first time that the total photoacoustic response in the presence of an impedance mismatch is a sum of the mismatch-free response and its time-delayed, time-reversed replicas, which may or may not overlap. We then show how to exploit the theoretical results for absorber reconstructions when the mismatch-free response does not overlap with its replicas (the case with short pulses) and for the more complicated situation when they do overlap (longer pulses or pulse-equivalent waveforms). These ideas are demonstrated via simulated absorber reconstructions, in particular with long, time-varying LFM pulses used with photoacoustic radars, as opposed to the short pulses typically used in conventional pulsed photoacoustics.

Theory

Definition of the 1D Boundary Value Problem

For a plane sample, of thickness L, immersed in a transparent fluid, the geometry of the problem is shown in Fig. 1.

Fig. 1

geometry of 3 layer problem.

In Fig. 1, c1, ρ1 are the speed of sound and density in the transparent fluid, respectively,c2, ρ2, μ, C are the speed of sound, density, optical absorption coefficient and specific heat of the absorbing layer. The back-propagated pressure response is denoted by p(z, t) in the first layer, corresponding to a measurement made in the reflection mode. The middle layer pressure response is denoted by p(z, t) and the transmission pressure response in the third layer is denoted by p(z, t). The governing equations for each layer are given in the frequency domain bywhere A(z) is a function of space that describes the geometry of the 1D absorber (depth dependence of the absorption coefficient). The corresponding boundary conditions are continuity of pressure and acceleration of the particles at both boundaries, which are given in the frequency domain respectively byIn this development, the non-unitary, angular frequency version of the Fourier transform is used.

Analysis of the Pressure Response

For the case of a finite layer, the pressure response must consist of only outwardly propagating waves in the first and third layers, but must have waves that propagate in both directions within the finite layer itself. Therefore, the pressures for the finite case are given byIn equation (4), k = ω/c for i = 1, 2 and B(z, k2) is the particular solution to the nonhomogeneous equation for (second equation in (2)), and thus depends on the precise form of . It is written as to clearly denote its dependence on z and k2 = ω/c2, rather than ω alone. It is noted that B(z, k2) does not have any dependence on c1, since c1 does not appear in the equation which defines it. For example, if the absorber is modeled for Beer's Law absorption so that A(z) = e−, then

The integration constants and , also depend on the angular frequency ω and are denoted as such. The functions e± are the outwardly propagating traveling wave solutions. The specific expressions for the integration constants are found by applying the boundary conditions at z = 0 and z = L (equations (3)). To simplify the expressions, the following definitions are used: first, the impedances are defined as Z = cρ and also the reflection and transmission coefficients are defined asNote that these satisfy T = 1 + R and that when there is no impedance mismatch between the two media (that is when Z1 = Z2), then R = 0, T = 1. Also note that since all impedances are positive numbers then , a fact that shall be useful later.

The unknowns and , can be found from applying the four boundary conditions in equation (3), which gives four equations in four unknowns. These can be solved and further simplified using the definition of reflection and transmission coefficients given in (6). For the back-propagating pressure-wave, z < 0, the pressure is given bywhereand G is the complex conjugate of G given byFor physical problems, the functions and therefore are real-valued.

Similarly, the middle-layer pressure is given byFinally, the transmission pressure isThe experimentally relevant back-propagation, , and transmission, , pressure waves will now be examined separately.

Analysis of the Back-propagation Pressure Wave

In this section, we demonstrate that under certain conditions, inversion algorithms that do not model an impedance mismatch may be used even if a strong mismatch is present. First, we note that the expression for can be written asHere, TF is the total transfer function of the system, which encapsulates the effects of the geometry of the absorber along with the impedance mismatch of the media via the reflection coefficient R and is given byThe transfer function components in equation (12) are, in turn, given byFrom the form of equation (11), it can be easily seen that interaction with the layer modifies the form of the input pulse via the total transfer function given in (12).

The total transfer function in equation (12) can be viewed as a combination of the two leading terms with exponentials, along with the R-weighted sum of the sub-transfer functions TF0 and TF1. The first term in the transfer function is simply a time delay while the second term is easily shown to be an increasingly attenuated sum of time-delayed reflections of the response within the absorber as seen by the power series expansion for given byHence, only the transfer function components TF0 and TF1 will change the shape of the input pulse since the effect of the other terms is only one of time-delays and attenuations. Equation (12) can be written as

For real functions, the effect of conjugation in the frequency domain is the same as time reversal in the time domain. The effect of the final term, in equation (16) is the sum of (i) a system-modified shape via the TF0 term and (ii) the negative, attenuated and time-delayed, time-and shape-reversed term via . Hence, equation (16) demonstrates that the total response is a sum of the mismatch free response and time-delayed, time-reversed versions of it. Simulations with equation (16) assuming Beer's law absorbers and various waveforms produce the classic N-shaped profile, [20], [24], now with the additional insight that the classic shape is the superposition of the TF0 and TF1 terms.

The duration of the TF0-controlled term depends on the duration of the input pulse and on the shape of the absorber via TF0 itself. The term is time reversed and clearly has the same duration, although it does not appear until time 2L/c2 later. Therefore, these two components may or may not overlap (and therefore sum) depending on the relative sizes of the duration of impulse response and the transit time through the absorber. Furthermore, TF0, the transfer function due to the shape of the absorber, would be the entire transfer function for the case of no-impedance mismatch (R = 0). Clearly, for small values of R, and TF0 alone can represent the effect of the absorber, neglecting the impedance mismatch. It is also clear that as R increases, the contribution of will not be negligible because the two terms will be of the same order of magnitude since .

When there is no temporal overlap between TF0 and , TF0 accounts for the compression wave and for the rarefaction wave. Therefore, because of the lack of temporal overlap, TF0 could be deduced from the compressive portion of the wave only–even for the case of a large reflection coefficient. This implies that for simple objects, neglecting the rarefaction portion of the wave gives the shape of the wave that would be obtained with the photoacoustic equation without the effect of the impedance mismatch (R = 0), even for the case of a large R. For the purposes of the inverse problem of reconstructing the shape of the absorber from measurements of the response, the compressive wave alone would yield sufficient information to reconstruct the shape of the absorber using the 1D photoacoustic equation with no impedance mismatch. This only follows because there is no temporal overlap between the TF0 and term and will occur whenever the duration of the response to TF0 is equal to or smaller than L/c2. For an impulse response, the duration of the response to TF0 should be at most L/c2, therefore this condition applies in the case of short-pulse photoacoustics. Once the duration of the response to TF0 is greater than L/c2, then the compressive wave does not represent TF0 alone. This may be the case with longer pulses and pulse-equivalent waveforms such as LFM chirps. These points are demonstrated via simulation in Appendix A.

Importantly, for the case of the delta-function input waveform, careful selection of the appropriate (compressive) portion of the measured response does permit the use of algorithms that assume no impedance mismatch - even if a strong impedance mismatch is present. While an actual delta function waveform is impossible to create, it would suffice if Δt << L/c2 for this case to hold. For longer duration waveforms, TF0 itself cannot be deduced from the shape of the compression or rarefaction wave. The implications of this for absorber spatial reconstruction algorithms is that for longer duration waveforms, it is imperative to take into account the impedance mismatch in such algorithms.

Analysis of the Transmitted Pressure Wave

In this section, the expression for pressure measurements made in the transmission mode is similarly analyzed to demonstrate that under certain conditions, inversion algorithms that do not model an impedance mismatch may be used even if a strong mismatch is present. It can be shown that the expression for can be written aswhereAs for the case of the back-propagated pressure wave, equation (18) demonstrates that the total response is a sum of the mismatch free response and time-delayed, time- and shape-reversed versions of it. In equation (18), there is no time delay between the and the R · TF0 term, although the term is a time-reversed version of the TF0 term (a consequence of complex conjugation in the frequency domain).

For an ideal delta-function impulse response ), the inverse Fourier transform will be the exact time reversal of and there will be no temporal overlap because both contributions start at the same time and proceed in different directions in time. In this case, selecting only the rarefaction or compression portion of the response would immediately yield the corresponding mismatch-free response, which could then be used in inversion algorithms that do not model the impedance mismatch - even if a strong impedance mismatch is present. This is a consequence of the delta-impulse nature of the input pulse which ensures non-overlap between the two transfer function contributions since they each proceed in different directions in time. While an actual delta function waveform is impossible to create, it would suffice if Δt << L/c2 for this case to hold.

However, for non-delta pulses of finite length the situation is different. As the input pulse widens from the ideal case of an input delta function , the two response terms and start to widen and overlap in time. The R-order term (proportional to TF0) starts at the earliest at and proceeds forward in time with a duration of at least L/c2 + Δt, where Δt is the duration of the input pulse. On the other hand, the 0th order term (proportional to ) starts at time and then proceeds ‘backwards’ in time for a duration of L/c2 + Δt to time . Both 0th and 1st order terms have the same duration but the overlap occurs because they do not start at the same time as was the case for the impulse response. This is the cause of the temporal overlap in the two component responses.

For long pulses or pulse-equivalent waveforms like LFM chirps, and a large R value, it would not be possible to isolate the mismatch-free response from a simple analysis of the measured response. Therefore, TF0 itself cannot be deduced from the shape of the compression or rarefaction wave. The implications of this for absorber spatial reconstruction algorithms is that it would be imperative to take into account the impedance mismatch in such algorithms.

Sample reconstructions using back-propagating waveforms

In this section, we demonstrate by way of example how the ideas developed in prior sections can be put in place for an absorber reconstruction from a simulated measured back-propagated pressure waveform. Equations (11), (12) give the exact theoretical relationship between the back-propagated pressure wave, input temporal waveform and absorber that give rise to the transfer function. We assume that measurements of the back-propagation pressure wave are time-gated to remove multiple reflections, which has the effect of gating out the second and higher terms of the series expansion of the denominator in equation (15) so that only the first term survives. The signal is then shifted to the origin in time to remove the effect of the time taken for the waveform to reach the detector from the sample and finally normalized (scaled) since the absolute amplitude of the signal contains no information. Then the equivalent relationship to equation (11) only contains the effect of the transfer function components and input waveform and can be represented bywhere we have used the superscript exp to denote a measured wave that has been processed to be ready for analysis (time-gated to remove reflections, shifted to the origin in time and scaled).

For the simple Beer's law absorber we considered above, for 0 ≤ z ≤ L, so is given by equation (5), is given by equation (9) and Eq. (13) for the transfer function component TF0 givesWriting the absorber function as and taking the spatial Fourier transform of so that the spatial variable z transforms to spatial frequency ω, givesWe note that the spatial Fourier transform of evaluated at ω = k2 givesIn other words, to within a proportionality constant, TF0 gives the spatial Fourier transform of the absorber. Using the fact that , and making use of the observation in (22), then for the simple Beer's law absorber, equation (19) can be written as

Using the ideas discussed in the previous sections, equation (23) is a simple way to arrive at the absorber profile from the measured back-propagated pressure wave. In the first instance, for a small reflection coefficient where R ≈ 0 equation (23) reduces toand the absorber profile can be obtained by inverse Fourier transformation of equation (24), which amounts to deconvolving the shape of the input pulse from the measured pressure wave. It is noted that the speed of sound in the absorber, c2, must be known a priori in order to obtain the size of the absorber, otherwise equation (24) will return only the shape of the absorber profile. In the special case that the input waveform is a delta function pulse or is of sufficiently short duration with respect to the transit time of the absorber to be approximated as a delta pulse then and equation (24) directly gives .

Reconstructions using the compressive portion of the wave profile

Now, to put into practice the ideas proposed in previous sections for a non-negligible reflection coefficient R, if the input pulse is a delta function pulse or is of short enough duration with respect to the transit time of the absorber to be approximated as a delta pulse, then most of the second term in equation (23) is contained in the rarefaction portion of the measured pressure wave. Therefore, taking only the compressive portion of the measured wave gives an equivalent expression to equation (24)–even for large values of R:In equation (25), the subscript R,C has been used to denote the compressive portion of the back-propagated pressure wave. In other words, even for the case of large values of reflection coefficient, under a suitably short-lived input waveform and by ignoring the rarefaction portion of the wave, the same reconstruction approach that can be used for a mismatch free case (equation (24)) can also be used for the case with a large impedance mismatch (equation (25)). Interestingly, Anastasio et al. noted that the deterioration in the quality of acoustic signals scaled with depth [23] and exploited this to improve reconstruction quality by using only the leading half of the collected acoustic signals. Their logic was that the trailing half of the signal corresponded to regions in the object that are farther than the object's center. The approach was heuristic since analytic inversion formulae for the reconstruction of images from half-time data functions had not been identified. They used expectation-maximization (EM) based reconstructions to demonstrate that improved reconstructions in acoustically heterogeneous media could be obtained only if the first half of the signal was used rather than the entire signal. What we have shown in equation (25) is that their heuristic observation indeed has a mathematical foundation.

On the same subject, Cox and Treeby [15] considered the effect of acoustic heterogeneity on photoacoustic reconstructions, focusing on a time-reversal approach. They showed that the boundary condition inherent in time reversal imaging can act as a reflector to trap the time reversed versions of waves scattered by acoustic heterogeneities, leading to a greater level of artifact in the image. They argued that the enforced time reversal boundary condition can “trap” artifacts in the final image, and that by truncating the data, or introducing a thresholded boundary condition, this artifact trapping can be mitigated to some extent. These observations were made specifically for time-reversal based reconstructions based on pulse (impulse excitation) photoacoustics. The ideas proposed in the present paper of time-gating to remove the multiple wave reflections and taking only the compressive portion of the pulse correspond to their suggestions of truncating and thresholding the data, with the difference that the comments made here are independent of the form of the reconstruction algorithm and do not rely on a delta-impulse input pulse.

However, it is important to note that once the input waveform is of sufficiently long duration that the contributions of TF0 and TF1 start to overlap and add so that they cannot be separated by simply removing the rarefaction tail of the wave, then use of equation (25) will produce errors. We illustrate these ideas by way of simulations.

For the simulations, we consider parameters to highlight the effect of large or small impedance mismatches. The critical parameter in our development is the reflection coefficient, R, which represents the degree of impedance mismatch. The actual material parameters themselves do not affect the shape of the wave profile other than indirectly via R. Therefore, material parameters were chosen so that the reflection coefficient would be small (close to zero, small impedance mismatch) or large (close to 1, large impedance mismatch) in order to investigate the effect of a small or large impedance mismatch on the results. For a water and steel interface, the parameters are μ = 200 m−1, L = 5 mm, c2 = 1500 m/s(water), c2 = 5790 m/s (steel), ρ1 = 1000 kg/m3 (water), ρ2 = 7480 kg/m3 (steel) which gives R = 0.993 and a transit time through the absorber of 0.86 μs. Steel was chosen for demonstrative purposes since the impedance mismatch with water is very large. For a water and blood interface, we replace the speed of sound and density of the absorber with c2 = 1570 m/s (blood) and ρ2 = 1060 kg/m3(blood) [25], which gives R = 0.0519 and a transit time through the absorber of 3.18 μs. This represents a typical case with a negligible impedance mismatch.

It is noted that biological tissue generally has much higher variability than water, blood and steel. For example, the speed of sound varies from approximately 600 m/s in lung tissue to circa 4,000 m/s in bone. In the brain, the speed of sound has been measured to about 1,550 m/s [26]. There are also large variabilities in the reported values of measurements of speed of sound for the same biological tissue [27]. Densities in biological tissue are also equally variable with ranges of 260 kg/m3 in lung tissue to 1120 kg/m3 in bone [28]. Again, there is variability in the reported measurements of densities for the same biological tissue. For example, reported values for bone differ from 1120 kg/m3 in [28] to 1850 kg/m3 in [29]. For the reported impedances of bone and water, the reflection coefficient would be R = - 0.66. For blood and lung tissue, the reflection coefficient would be R = 0.83. Note that metals are often used to pin bones after a break in the bone.

As discussed above, we first consider the case of a blood-water interface so that R = 0.0519, representing a typical case with a negligible impedance mismatch. For a delta function impulse, in the simulations and we additionally consider a short rectangular pulse of duration 1.67 μs. The 1.67-μs value was chosen in order to enable comparison with other simulations used for experimental validation in this paper (Appendix A). Since the transit time for this absorber is 3.18 μs, this short pulse is not short enough to be considered an effective delta pulse for this absorber. However, R is small enough to be very close to negligible and so it is anticipated that there would not be major error in ignoring the impedance mismatch. Fig. 2 uses equation (19) to model the response for both the delta pulse and the short pulse for both cases where the full transfer function is taken and for the case where R = 0 is set in the simulations. Since Fig. 2 (a) is the impulse response, it is expected that its temporal shape is the same shape as the spatial profile of the absorber given by , and this is what is shown in the simulations: a Beer's law absorber of length 5 mm amounts to an exponentially decaying temporal profile of length L/c which is the transit time of 3.18 μs. It is noted that, as would be expected, the longer pulse broadens and changes the shape of the response. Also as was expected, the simulations where R = 0 is considered in the total transfer function do not differ appreciably from the full transfer function case.

Fig. 2

Back-propagation pressure responses for blood-water for (a) delta pulse and (b) short (1.67 μs) pulses. Parameters are μa = 200 m−1,c1 = 1500 m/s (water), c2 = 1570 m/s (blood), z = -0.02 m, Δt = 1.67 μs, ρ1 = 1000 kg/m3(water), ρ2 = 1060 kg/m3 (blood).

The absorber is then reconstructed by taking the inverse Fourier transform of equation (25) using the simulated forward response for experimental pressure data. The compressive portion of the response is found by taking the positive portion of the response (in this case, not appreciably different from the total response). The reconstructed absorber profiles are shown in Fig. 3. It can be seen that in both cases, the absorber profile is successfully reconstructed by inverse Fourier transformation of equation (25), even for the case where the pulse is not a delta impulse. It is also noted that the case with a short pulse produces more numerical artifacts that the delta pulse case. This is a result of the deconvolution process, which was numerically implemented as division by and thus tends to amplify high frequency effects. This is a numerical artifact of the approach taken and an appropriate low-pass filter design for the deconvolution process would likely produce better results than simple division by .

Fig. 3

Absorber reconstructions using compressive part of delta pulse and short (1.67 μs) pulse.

The same approach can now be considered in the more interesting case of a water-steel interface where the reflection coefficient is a non-negligible R = 0.993, representing a typical case with a large impedance mismatch. The cases considered are those of a delta pulse , a short square pulse (duration 1.67 μs) and a shorter square pulse (duration 0.1 μs). The 1.67-μs value was chosen in order to enable comparison with other simulations in this paper (Appendix A) and 0.1 μs was chosen since it is one order of magnitude smaller (less than a tenth) of the value of the 1.67 μs value and thus represents a comparatively short pulse. As before, the spatial profile of the absorber given by , with a chosen absorber of length 5 mm, implying an exponentially decaying temporal profile of duration L/c. The transit time, L/c, for this absorber is 0.86 μs so that the short 1.67 μs pulse is not short enough to be considered a delta function pulse.

Fig. 4 shows the response for the delta pulse and the short (0.1 μs) pulse, and Fig. 5 shows the response for the delta pulse and the longer (1.67 μs) pulse. Both figures show the results for both the case where the full transfer function is taken into account and for the case where R = 0 is set in the simulations. From Fig. 4, it can be seen that the shorter 0.1 μs pulse is indeed short enough to emulate a delta pulse and furthermore that the compressive portion of these ‘short enough’ pulses matches the R = 0 response.

Fig. 4

Back-propagation pressure responses for a water-steel interface, (a) response to delta and (b) response to shorter (0.1 μs) pulse. Parameters for water and steel are μa = 200 m−1,c1 = 1500 m/s (water), c2 = 5790 m/s (steel), z = -0.02 m, Δt = 1.67 μs, ρ1 = 1000 kg/m3(water), ρ2 = 7480 kg/m3(steel).

Fig. 5

Back-propagation pressure responses for a water-steel interface, (a) response to delta and (b) response to shorter (1.67 μs) pulse.

However, it can be seen from Fig. 5 (b) that the 1.67 μs pulse shows a significant broadening of the compressive peak and thus a deviation from the impulse response. More importantly, the compressive portion of the 1.67 μs pulse response deviates significantly from the R = 0 response, which implies that significant errors can be expected if this compressive portion of the wave is used in any reconstructions. The reconstructions will confirm this hypothesis. A practical conclusion from the simulations shown in Fig. 4 and Fig. 5 is that broadened pulses and pulse-equivalent waveforms can behave like ideal delta functions. It would be expected that pulse-equivalent waveforms that have a bandwidth similar to a ‘short enough’ pulse will yield photoacoustic signals equivalent to a short laser pulse. For longer pulses (equivalently, pulses with a smaller bandwidth), excessive broadening results in the case where there is an overlap between outgoing and back-reflected components of the signal and the full wave profile must be used in the reconstruction, as will be shown later.

Fig. 6 shows reconstructions performed by using the compressive portion of the very short 0.1 μs pulse (Fig. 4) in equation (25). The results are similar to those for the small R case: the delta impulse response produces a very good reconstruction and the reconstruction with the shorter pulse still produces a reasonable result although with a few more oscillations that can be attributed to the process of deconvolving the input pulse out of the results, due to the finite pulse bandwidth. This confirms the theory that was presented in the first part of the paper–that for impulses or ‘short enough’ pulses, the compressive portion of the measured response can be used in the same reconstruction formulas, regardless of the size of the impedance mismatch.

Fig. 6

Absorber reconstructions using compressive part of delta pulse and very short (0.1 μs) pulse.

However, when the pulse is no longer short enough to be considered an impulse, then the same approach fails. Reconstructions using the same approach of taking the compressive portion of the measured response to a short (1.67 μs) pulse are shown in Fig. 7. Clearly, the results are not meaningful and this confirms the prior statement that serious errors should be expected when using the mismatch-free reconstruction formulas with non-small pulse responses of absorbers with large impedance coefficients

Fig. 7

Absorber reconstructions using compressive part of delta pulse and short (1.67 μs) pulse.

In conclusion, with large impedance mismatches and an insufficiently short pulse, the full wave profile consisting of compressive and rarefactive portions must be used in the reconstruction. This ability to use only the compressive portion of the wave in reconstructions could be viewed as one of the advantages of pulsed photoacoustic imaging since typical durations of pulses are in the nanosecond range and therefore should be sufficiently short to ensure no overlap between the compressive and rarefactive portions of the pulse. However, ultrasonic transducers with bandwidth sufficient to faithfully record the true shape of the pulse would also be required.

Reconstructions using the full wave profile

Returning to equation (23), we denote so equation (23) givesTherefore summing and rearranging givesEquation (27) is used for reconstructions, utilizing the full simulated back-propagated pressure wave response and the short 1.67 μs pulse (same as used in Fig. 7). Results are shown in Fig. 8. Comparison of Fig. 7 and Fig. 8 makes it clear that for a pulse that is insufficiently short to be considered a delta function, the full pressure response must be used in the reconstruction in order to obtain acceptable results. More importantly, the reconstruction algorithm itself must be modified to incorporate the effects of the impedance mismatch.

Fig. 8

Absorber reconstructions using full response to short (1.67 μs) pulse.

Application to LFM waveforms

These ideas can also be adapted to (longer-pulse) shaped waveforms. In particular, recent work has shown that chirped modulation waveforms achieve dramatic signal-to-noise increases in the photoacoustic measurements [30]–[32]. Using the same parameters as in the previous steel-water example, the 1.67 μs pulse is now modulated so thatEquation (28) describes a linearly modulated pulse that turns on and off at times t0 and t respectively and consists of a chirp with starting angular frequency ω0 and a chirp rate of β. Fig. 9 shows simulations using the same parameters as in the previous example with a short 1.67 μs-duration chirp modulated according to (28) where ω0 = 106rad/s and β = 100 · 106rad/s2. The high frequency modulations can be seen in the response of Fig. 9(b) for both the full and mismatch-free case.

Fig. 9

Simulations for a water-steel interface, (a) Linearly modulated input chirp (b) Back-propagated photoacoustic response to (1.67 μs) short chirped pulse.

Following in the same manner as the previous example, reconstructions can be performed by taking only the compressive portion of the response and using equation (25). These results are shown in Fig. 10. Clearly, the now-linearly modulated pulse is still insufficiently short to be considered a delta function, and the full pressure response must be used in the reconstruction in order to obtain acceptable results. Equation (27) is the full-model response and can be used for optimized reconstructions. The absorber is again well reconstructed in this case, leading to reconstructions very similar to Fig. 8.

Fig. 10

Absorber reconstructions using the compressive part of delta pulse and short (1.67 μs) linearly chirped pulse.

Telenkov and Mandelis discuss the increased signal to noise ratio that could arise when the Photoacoustic Radar (PAR) is used. In PAR, the output pressure pulse is cross-correlated with the input temporal waveform [30] to obtain the PAR waveform. We show here that the same procedure for reconstructing the absorber profile using the full model can also be applied to the PAR case. Using the fact that cross-correlation in the time domain implies multiplication by the complex conjugate in the frequency domain, returning to equation (23), it follows thatSuppose that now is defined as the ratio of the cross-correlation of pressure with input waveform (PAR pressure response) to the auto-correlation of that same waveform:Using (30), it then follows that equation (27) still holds with this modified definition of . Using the same steel-water model with the chirp of Eq. (28) and the same parameters used to obtain Fig. 9, the absorber can be reconstructed using equation (27) and the PAR pressure-response. This approach gives an excellent reconstruction, similar to that shown in Fig. 8.

In summary, when using the full pressure pulse response in reconstructions with the approach outlined in this paper, the photoacoustic radar pressure response can be successfully used in absorber reconstructions in the place of the pressure response itself. This is particularly useful since, as discussed in [31], [32], a major goal of PAR is to achieve a photoacoustic pressure response with higher signal-to-noise ratio than the pressure pulse response alone.

Conclusions

In this paper, the implications for the photoacoustic effect of an impedance mismatch between an optical absorber and its surroundings were examined with a view to obtain physical insights into circumstances that permit inversion algorithms to be adapted in cases where impedance mismatches are large and/or input pulses are not short. This was done via an understanding of how the absorber information carried by the components of the generated pressure wave can be used for signal inversion and absorber depth and shape reconstruction. It was demonstrated that the total transfer function in the presence of an impedance mismatch (without taking transducer transfer function-generated broadening into account) can be constructed from the mismatch-free case. Importantly, it was shown that under conditions of a delta-function input waveform (or one of small enough duration), although the response of the impedance mismatched system differs non-negligibly from the mismatch-free case, a careful analysis of the measured response can immediately yield the corresponding mismatch-free response, which could then be used in inversion algorithms that do not model the impedance mismatch. The implication of this is that any existing reconstruction algorithm used for the mismatch-free case can be adapted for use in the presence of a strong impedance mismatch by appropriate conditioning of the measured pressure response. These ideas were demonstrated via some simple reconstructions using simulated back-propagating pressure responses. Once the duration of the input waveform is large enough to be of the same order of magnitude as the transit time in the absorber, the preceding approach cannot be taken and the system response must be modeled as a combination of both 0th and 1st order reflection coefficient terms. This implies that any inversion algorithm must be adapted for the presence of the impedance mismatch and the full-wave pressure response must be used in the reconstructions. Mismatch-dependent reconstructions using simulated back-propagated pressure responses were performed to illustrate these ideas. A key finding from the simulations is that non-pulse waveforms should always consider a full impedance mismatch analysis due to the possibility of temporal overlap between purely compressive and rarefactive signal components. Specifically, the analysis in this paper addresses LFM waveforms that are not true delta functions but produce photoacoustic signals that are similar to those produced by a true delta pulse, especially within the limits of the measurement instrumentation. Reconstructions were successfully performed with chirped and cross-correlated input waveforms, demonstrating that these ideas could be used in combination with signal processing techniques of photoacoustic radars.

Conflict of Interest statement

The authors declare that they have no conflicting interests.