Development of a heterodyne speckle imager to measure 3 degrees of vibrational freedom

A laboratory system has demonstrated the measurement of three degrees of vibrational freedom simultaneously through heterodyne speckle imaging. The random interference pattern generated by the illumination of a rough surface with coherent light can be exploited to extract information about the surface motion. The optical speckle pattern is heterodyne mixed with a coherent reference. The recorded optical data is then processed to extract three dimensions of surface motion. Axial velocity is measured by demodulating the received time-varying intensity of high amplitude pixels. Tilt, a gradient of surface displacement, is calculated by measuring speckle translation following extraction of the speckle pattern from the mixed signal. This paper discusses the laboratory sensor concept, signal processing, and experimental results compared with numeric simulations. © 2016 Optical Society of America OCIS codes: (120.0120) Instrumentation, measurement, and metrology; (280.0280) Remote sensing and sensors. References and links 1. R. Jones and C. Wykes, Holographic and Speckle Interferometry, Vol. 6 (Cambridge University, 1989). 2. P. Hariharan, Basics of Interferometry (Academic, 2010). 3. P. Jacquot and J. M. Fournier, Interferometry in Speckle Light: Theory and Applications (Springer Science and Business Media, 2012). 4. M. L. Prado, N. Bolognini, and M. Tebaldi, “Analysis of tilt by modulated speckles generated with a double aperture pupil mask,” Opt. Commun. 338, 473–478 (2015). 5. B. Redding, A. Davis, C. Kirkendall, and A. Dandridge, “Measuring vibrational motion in the presence of speckle using off-axis holography,” Appl. Opt. 55(6), 1406–1411 (2016). 6. T. W. Du Bosq and E. Repasi, “Detector integration time dependent atmospheric turbulence imaging simulation,” Proc. SPIE 9452, 94520B (2015). 7. P. Reu, D. P. Rohe, and L. D. Jacobs, “Comparison of DIC and LDV for practical vibration and modal measurements,” Mech. Syst. Signal Pr. (in press). 8. B. Weekes and D. Ewins, “Multi-frequency, 3D ODS measurement by continuous scan laser Doppler vibrometry,” Mech. Syst. Signal Pr. 58, 325–339 (2015). 9. J. M. Reyes and P. Avitabile, “Use of 3D scanning laser vibrometer for full field strain measurements,” in Experimental Techniques, Rotating Machinery, and Acoustics, Vol. 8 (Springer International, 2015), pp. 197–209. 10. H. Yan, H. Z. Duan, L. T Li, Y. R. Liang, J. Luo, and H. C. Yeh, “A dual-heterodyne laser interferometer for simultaneous measurement of linear and angular displacements,” Rev. Sci. Instrum. 86, 123102 (2015). 11. A. Pouya and M. Alireza, “Laser Doppler interferometry (LDI) to obtain full stiffness tensor: a case study on a deformation zone in Sweden,” in ASEG Extended Abstracts, (2015), pp. 1–4. 12. S. Pasquet, L. Bodet, Q. Vitale, F. Rejiba, R. Gurin, R. Mourgues, and V. Tournat, “Laser-Doppler acoustic probing of granular media with varying water levels,” Phys. Procedia 70, 799–802 (2015). 13. J. G. Chen, R. W. Haupt, and O. Buyukozturk, “Operational and defect parameters concerning the acoustic-laser vibrometry method for FRP-reinforced concrete,” NDTE Int. 71, 43–53 (2015). #258856 Received 8 Feb 2016; revised 25 Mar 2016; accepted 4 Apr 2016; published 8 Apr 2016 © 2016 OSA 18 Apr 2016 | Vol. 24, No. 8 | DOI:10.1364/OE.24.008253 | OPTICS EXPRESS 8253 14. J. Bencteux, “Holographic laser Doppler imaging of pulsatile blood flow,” J. Biomed. Opt. 20, 066006 (2015). 15. F. Tenner, D. Elz, Z. Zalevsky, and M. Schmidt, “Optical tremor analysis with the speckle imaging technique,” J. Imaging Sci. Techn. 59, 10402 (2015). 16. X. Lai and H. Torp, “Interpolation methods for time-delay estimation using crosscorrelation method for blood velocity measurement,” IEEE T. Ultrason. Ferr. 46(2), 277–290 (1999). 17. J. W. Goodman, Introduction to Fourier Optics (Roberts and Company, 2005). 18. H. J. Tiziani, “Analysis of mechanical oscillations by speckling,” Appl. Opt. 11(12), 2911–2917 (1972). 19. K. ODonnell and E. Mendez, “Experimental study of scattering from characterized random surfaces,” J. Opt. Soc. Am. A 4(7), 1194–1205 (1987). 20. L. F. Rojas-Ochoa, D. Lacoste, R. Lenke, P. Schurtenberger, and F. Scheffold, “Depolarization of backscattered linearly polarized light,” J. Opt. Soc. Am. A 21(9), 1799–1804 (2004). 21. J. W. Goodman, Speckle Phenomena in Optics: Theory and Applications (Roberts and Company, 2007). 22. B. K. Park, O. Boric-Lubecke, and V. M. Lubecke, “Arctangent demodulation with dc offset compensation in quadrature doppler radar receiver systems,” IEEE T. Microw. Theory 55(5), 1073–1079 (2007). 23. G. Sendra, H. Rabal, M. Trivi, and R. Arizaga, “Numerical model for simulation of dynamic speckle reference patterns,” Opt. Commun. 282(18), 3693–3700 (2009). 24. J. C. Dainty, Laser Speckle and Related Phenomena, Topics in Applied Physics, Vol. 9 (Springer-Verlag, 1975).


Introduction
The independent observation of Doppler phase for velocity and speckle shifts for tilt are well documented in the literature [1][2][3][4].Heterodyne Doppler vibrometry measures object velocity in the axial direction of the interrogation beam.Speckle image correlation measures displacement gradients, i.e., out of plane tilts, for two rotational dimensions, Fig. 1.This paper presents a model and laboratory system that combines these techniques, allowing simultaneous observation of three degrees of freedom using a single interrogation beam.It uses a heterodyne process to gather Doppler phase information where a coherent reference signal is modulated by a local oscillator and mixed with a measurement leg.The mix signal contains the Doppler phase, but because a focal plane array is used, spatial redundancy mitigates speckle dropout in laser Doppler sensing [5].Additionally, the focal plane array captures the speckle pattern, which is used to extract information about object tilting.Heterodyne gain reduces laser power requirements for speckle imaging, while faster sampling minimizes the impact of atmospheric turbulence [6].Prior work has quantified multiple degrees of freedom using multiple single beam vibration sensors [7][8][9][10].The ability to use a single beam allows for applicability to a number of geophysical [11,12], industrial [7,13], and medical [14,15] applications.The model below shows that information related to tilt and Doppler is present at an image plane; this information can be used to estimate axial velocity and tilt.The axial velocity is extracted using in-phase and quadrature digital mixing.All pixels contain axial velocity from a single object location.The axial velocity information is in theory equal on all pixels, but experimental results indicate noise reductions are possible by averaging pixels.
The tilt measurement begins with the heterodyne signal as well, but unlike axial velocity, multiple pixels are required to track speckle shift that is proportional to tilt.In our technique, the spatial speckle pattern is recovered from the heterodyne mix signal by computing the temporal envelope of each pixel.The speckle shift is found by cross-correlating sequential speckle images and computing a sub-pixel shift amount from the cross-correlation matrix [16].

Theory
The following describes an analytic representation of a heterodyne speckle imager and describes how motion in three degrees of freedom alter amplitude and phase of a speckle image.The spatially varying expressions follow previous work with some modification [17,18].The model calculates the electric field at an image plane following illumination of a rough surface.Spatial information at a detection plane is calculated using the Fraunhofer diffraction equations.These calculations are typically independent of time, but in this case motion is introduced as a phase shift due to path length change near the source.The optical frequency and a local oscillator are added later as phase shifts, as these phases do not affect the Fourier transform in the Fraunhofer equations.
In this model, an object plane is described in a {ξ , η, z} coordinate system, shown by Fig. 2, where ξ and η are in the plane and z is orthogonal to the plane.The model accounts for two types of motion: out of plane translation in the axial direction z and tilting of the plane in two directions.Out of plane motion results in a Doppler phase shift and tilting results in a translation of the speckle image [19].A single dimension of tilting is presented to simplify illustration, however, a second tilting degree of freedom is added to the final results.

Image Plane
Object Plane Unperturbed The moving object has a pivot point that is allowed to translate in the z direction, but is otherwise located at a fixed position ξ = ξ o .The position of a scattering point at {ξ P , 0, 0}, originally on the z = 0 plane, will be moved to a new location.This motion extends the path length of reflected light; the path length is used to determine the speckle image resulting from these two locations in a method paralleling Tizianis [18].
The variation in path length ∆P due to tilt and translation is a combination of the incoming and outgoing light.The incoming light is assumed to be coherent and planar at an angle β relative to the object plane.The image plane is assumed to be sufficiently far such that outgoing light is orthogonal to the object plane.The path variation is geometrically determined, (1) where θ is angular change due to the tilting surface to which a small angle assumption is applied: sin (θ ) ≈ θ , cos (θ ) ≈ 1 − θ 2 /2, and θ 2 θ .This assumption leads to a simplified path length, In the description below, the direction of illumination β = 0 and hence the Doppler velocity is assumed to be aligned with the z axis.The path length reduces to The electric field at the image plane before tilting U 1 (x, y) is described by the Fraunhofer diffraction equation, where u o (ξ , η) is the complex electric field amplitude resulting from the scattering object, λ is the optical wavelength, k = 2π/λ is the wave number, z is the distance to image plane that is assumed to be in the far field.For simplicity, the propagation coefficients outside the integral will be defined as κ = e jkz jλ z e jk 2z (x 2 +y 2 ) .
Substituting the variables we rewrite equation ( 5) in a compact form.
where F [ f (ξ , η)] (x, y) is the Fourier integral operator.Once the object is translated and tilted into a new position, an additional phase shift ∆P is added that represents the optical path length change The electric field at the image plane resulting from the tilted surface can now be described by Solving the Fourier transform yields where M ( f x , f y ) is the Fourier transform of the complex electric field of the scattering object u o (ξ , η).The above equations describe the field for two fixed positions.It is informative to consider the expression of the tilted image field U 2 (x, y) as a phase modified version of the original field without tilting U 1 (x, y), The tilted image field may now be expressed as a function of the undisturbed image field.
This shows that the tilted field variable is a spatially shifted and phase modified version of the original field.The irradiance of this field is a shifted version of the original.This derivation provides an analytic model for the amplitude and phase of light scattered from a tilted object.If the optical strength of the light is sufficient for measurement, this shift may be calculated and used to estimate the tilt of the surface θ .The objective of the work is to simultaneously measure velocity of the surface along the z-axis by incorporating a Doppler vibrometer technique.Thus, using a traditional heterodyne process to extract the phase of light, while still measuring the optical shift associated with object tilt.
The heterodyne process of extracting Doppler phase shifts and field shifts due to tilting is described below.Time will be added to the field variable expressions to capture changes to the phase of light.The angle of tilt θ becomes a function of time θ (t), that is associated with the moving object.A second, orthogonal dimension of tilt φ (t) is also added.The field variable described in equation (11) will be referred to as the measurement and a subscript "m" will be used to distinguish it from the heterodyne reference beam with a subscript "r."The image field is where ω o is the optical frequency, η o is the pivot in the second dimension and The term M (x (t) , y (t)) is dependent on time because the angles θ (t) and φ (t) are slowly changing.Since M (x (t) , y (t)) is complex, it can equally be represented as a magnitude |M| and phase α M .In this form the optical mixing equations will be easier to simplify.
A reference field U r is coherent with the measurement leg and shifted in frequency by ω LO , The reference field uses the same propagation distance z as the measurement field to match the radius of curvature, but is not subject to phase shifts caused by the object's motion.The total field U t at the image plane is the superposition of the measurement and reference fields The irradiance on the detector is The irradiance is the sum of three terms: 1) irradiance of the reference R 2 , 2) irradiance of the measurement |M| 2 and 3) the mix term 2R|M| at the local oscillator frequency.The mix term magnitude is related to the magnitude of the speckle pattern and therefore shifts spatially as this pattern shifts.R is assumed to be greater than |M| such that R|M| will be larger than |M| 2 and is used to determine the speckle shift as described in the next section.The magnitude of the measurement field is not affected by the axial motion in the z direction, but would be affected by translation orthogonal to the incident beam.
The local oscillator frequency ω LO is phase modulated.The phase is The term d 4 is the component of the Doppler shift due to translation of the object in the z dimension.Other terms in this phase are related to crosstalk between the axial and tilting motion.The terms containing θ and φ are Doppler shifts due to the tilting motion.It occurs in this expression because the pivots do not intersect the optical axis and the illumination is assumed to be centered at the origin.In other words, this term represents the component of motion in the z direction due to tilting, as applied over the distance between the pivot point and the center of illumination.
The terms described in the previous paragraph do not depend on image plane position.In contrast, α M is a function of x and y.Recall that α M is the phase of the measurement speckle pattern and shifts in a manner similar to the magnitude of the speckle pattern.This shift causes individual pixels to have a small phase modulation, thus adding noise to a Doppler measurement.

Laboratory concept
The laboratory sensor, represented below in Fig. 3, is similar to a conventional heterodyne Doppler vibrometer with the exception of optics for defocusing the measurement path and modification of speckle size.In this case, the focal plane array is used to measure spatial irradiance patterns instead of a single detector element.The laser source is a Thorlabs 20mW, 632nm linearly polarized HeNe.A 90:10 splitter generates a measurement and reference leg.The measurement and reference legs are frequency shifted with acousto-optic modulators (AOM) by 40 MHz and 40.01 MHz, respectively.After mixing this modulation generates a 10 kHz carrier frequency.
The measurement beam is focused onto a diffuse target that is mounted on a scanning mirror, Newport FSM-200, of which the reflective surface property is not used.The scanner is driven with an oscillating voltage to provide a dynamic tilt in one dimension only.
The beam is incident on the diffuse target with a small angle of incidence β ≈ 0 • .The return light is collected by a 50mm diameter, 500mm focal length (f/10) bi-convex lens that is aligned with the illumination location and oriented parallel to the diffuse target.A second lens, 25mm diameter and 250mm focal length (f/10), projects the speckle pattern onto the focal plane array while providing adjustment to the degree of defocus.The measurement beam becomes partially depolarized [20,21], following scattering from the target which interferes incoherently on the detector.To remove this, both legs pass through polarizing filters before alignment at the 50:50 splitter so that uniformly polarized light is mixed at the detector.In addition, the reference beam is expanded, filtered, and collimated to create a Gaussian wavefront.Only the full width at half maximum (FWHM) diameter mixes with the speckle pattern on the array.The camera, Vision Systems Phantom v341, uses a focal plane array that records in a color format at 128x128 pixels.Only the red pixels are processed, creating an effective 64x64 pixels at 30k frames/s.
The spatial, temporal data cube produced for each measurement is processed to calculate tilt in the x and y directions as well as axial velocity in the z direction.For the calculation of tilt, it is necessary to estimate the measurement leg's speckle pattern amplitude from the measured irradiance.This is obtained from the mixed components by temporally high-pass filtering to keep only the contents centered at the local oscillator frequency ω LO . #258856 Filtering removes the R 2 and |M| 2 terms from equation (21).Since the spatial variations of the measurement |M| are of interest, the instantaneous-amplitude envelope 2R|M| is estimated with a Hilbert transform H.

2R|M(x, y,t)|
The process removes removes instantaneous-phase fluctuations associated with the 10kHz carrier.
In practice, R is constant and larger than |M|; it provides heterodyne gain and allows for detection of the measurement field in cases where the intensity |M| 2 is too small to detect.The resulting speckle pattern at each frame is cross-correlated with a prior time frame to determine the number of pixels shifted.The object's tilt difference ∆θ is proportional to the spatial shift s between two subsequent frames, Vibrating objects generally have small amplitudes of tilt, this in turn produces image shifts that are a fraction of a pixel.To accurately determine these sub-pixel shifts, parabolic interpolation of the cross-correlation function is used to estimate the image shift.This method is applied in two dimensions, but results are shown only for a single tilt dimension under excitation.In the setup we use a random scattering surface that produces a stochastic speckle pattern.As a result, the cross-correlation technique is not exact when estimating the peak and corresponding spatial shift.
The calculation of velocity requires demodulation of the Doppler shifted signal relative to the local oscillator at individual pixels.The demodulation technique uses a standard arctangent demodulation [22] of digitally produced in-phase and quadrature signals.To generate I demod and Q demod , the measured irradiance is high-pass filtered to remove the R 2 and |M| 2 terms from equation (21).The phase information carried at 10kHz is shifted to base-band through digital mixing and filtered to remove frequencies aliased down from the 2ω LO component.
The in-phase and quadrature terms are then used to estimate the Doppler phase The velocity is calculated using the Doppler phase The vibration of interest is corrupted by speckle-phase noise ∂ α M /∂t when the tilt magnitude is large, but in our experiments this was not the case.Removal of the speckle-phase noise may be possible given that the speckle shift amount s is known.
Since the optical setup is defocused, each pixel should have the same z direction velocity.This redundant information is used to reduce noise by averaging across pixels.We are able to further improve the velocity estimation by selecting pixels that have high optical strength.This is done by revisiting the reconstructed speckle pattern and selecting regions where the measurement irradiance is high; a technique that mitigates speckle fadeout noise normally present on single pixel detectors.Similarly when calculating speckle shifts, 30k frames/s is oversampled relative to the vibration frequencies used in our experimentation; therefor, filtering and downsampling in time reduces noise.

Experimental results
Initial experiments were performed to validate the sensor concept and measure linearity.Ground truth was measured using a Polytec PSV-400 scanning laser Doppler vibrometer, which provided vibration amplitude and phase maps of the target surface to confirm that a linear velocity profile and uniform tilt were present.This also provided the exact location of the hinge point of the diffuse target where axial velocity is zero.A single axis of the target was excited by four simultaneous frequencies: 159, 190, 254, and 318 Hz. Figure 4 shows the z direction velocity amplitude of the target at 159 Hz, as measured by the Polytec and laboratory sensors.The laboratory sensor accurately determines axial velocity.In addition, the velocity profile shows the measurement surface is linear, which ensures a uniform tilt across the surface.Figure 5 shows the z direction velocity measurements for the laboratory sensor for each frequency (5Hz bandwidth) and spatial position.The z direction velocity is related to displacement and the tilt is this displacement as a function of position.However, for direct comparison with z direction velocity, the tilt measurements are expressed as a tilt velocity.The slope of each line estimates the tilt velocity of the surface for each frequency.The diffuse target hinge is located at 13mm.At this position, velocity goes to zero, however tilt velocity remains consistent along the surface, shown by the corresponding tilt velocity in Fig. 6.The tilt is proportional to pixels shifted, calculated between sequential frames.Here we divide by the per frame period to arrive at pixels per second, which is proportional to axial velocity for this object.The scaling to tilt velocity requires an estimation of the magnification factor of the optical system which has not been determined at the time of writing.In Fig. 6, an amplitude dependent variance is apparent on the pixel shift calculation.The velocity measurements indicate the surface is uniformly tilting, i.e., measurement of tilt at each spatial location should be identical.The variance is likely due to the accuracy of the crosscorrelation and parabolic interpolation routine.The routine may not be ideal for speckle images, in which new speckle patterns are generated from each spatial position on the target that is a random scatterer.A comparison with a numeric model, below, shows a similar variance.

Explanation of tilt variance
A numeric model is used to investigate stochastic variance in the calculated pixel shift using cross-correlation.It is hypothesized that a change in speckle pattern due to observation of a new surface results in variance on the tilt calculated.
This model, based on equation ( 8), generates an electric field at the image plane following illumination of a dynamic rough surface.Spatial information is calculated using discreet Fourier transforms, dependent on Fresnel diffraction equations [23].This transform is independent of time, but motion is introduced as a dynamic phase shift due to path length change and applied to the object field, u o (ξ , η), before the transform.The optical frequency is removed and the analysis is done without a heterodyne signal.Figure 7 shows speckle patterns acquired experimentally and generated with the numeric model.Multiple random instantiations of the rough surface were generated.The surface was then modulated with identical frequencies and velocities used in the experimentation.Figure 8 compares the calculation of image translation for modeled and experimental data.The experimental mean and variance were generated from the data shown in Fig. 6.Qualitatively, the variance on the pixel shift calculated using modeled data agrees with the experimental measurements, indicating a primary source of variance is due to the ability of the cross-correlation routine to accurately resolve sub-pixel shifts from stochastic patterns.Initial work with the numeric model shows that the tilt calculation is sensitive to the speckle size.For an objective speckle pattern, the speckle size is a function of the illumination area, wavelength, and distance to the observation region, represented by where A s is the correlation area and A o is the area of illumination.The square root of A s provides a reasonable approximation of speckle width [24].In the model, the incident beam size was adjusted to modify the average speckle size.The complex electric field from Eq. 5, u o can be replaced by where ζ is a spatial random phase operator and A o is a binary image mask.Figure 9 shows the statistics of pixel shift for 10 iterations of each speckle size.As speckle size decreases, the accuracy of the cross-correlation peak fitting routine improves, reducing variance.The experimental data was acquired using a comparative speckle size of 60µm, indicating it is possible to improve the estimation of tilt by reducing speckle size.

Conclusions
A heterodyne interferometric sensor concept has been generated which allows for the simultaneous observation of motion in three degrees of freedom.Axial velocity can be determined in the z direction and tilt can be determined about the object's x axis and y axis.Analysis shows the components of motion manifest as temporal (phase) or spatial variations when a heterodyne signal is measured.At present each component is extracted individually.Axial velocity produces a phase shift carried on the local oscillator and can be readily extracted using IQ demodulation.Surface tilt generates spatial speckle-pattern shifts in the measurement field and the heterodyne mix irradiance.Tilt is proportional to sub-pixel shifts in the speckle pattern.
Experimental irradiance data collected from a tilting, diffuse scattering target prove the technique is capable of measuring the motion components; the velocity demodulation matches calibrated sensors and tilt velocity calculations qualitatively agree with numeric models.The measurements and model indicate variance in the tilt estimation due to observation of new random surfaces and accuracy of the cross-correlation peak fitting routine.The variance may be reduced by optimizing the speckle size, as shown by the numeric model.
The observation of multiple degrees of freedom allows for a more complete view of arbitrary surface motion in situations where the surface has independent, axial and rotational motion components.The laboratory sensor concept provides an initial approach for quantifying these components simultaneously using a single beam.

Fig. 1 .
Fig. 1.Components of motion observed by heterodyne speckle imager.Velocity is measured in the axial direction z of the interrogation beam.Tilting motion is measured about the x and y axes.

Fig. 4 .
Fig. 4. Dynamic diffuse target velocity profile measured by a commercial Polytec LDV and the laboratory sensor.A single frequency, 159Hz, is shown.

Fig. 5 .
Fig. 5. Axial velocity amplitude at various positions on the diffuse target.Four simultaneous frequencies were present on the dynamic measurement surface.

Fig. 6 .
Fig.6.Pixel shift velocity corresponding to tilt velocity, acquired simultaneously with axial velocity.Four simultaneous frequencies were present on the dynamic measurement surface.