Estimating sound intensity from acoustic data captured byparallel phase-shifting interferometry
Fumihiko IMAEDA(1) , Kenji ISHIKAWA(2)∗ , Kohei YATABE(3) , Yasuhiro OIKAWA(4)(1)

Department of Intermedia Art and Science, Waseda University, Japan, emergency-utsuho@asagi.waseda.jp(2)

Department of Intermedia Art and Science, Waseda University, Japan, k-ishikawa@fuji.waseda.jp(3)

Department of Intermedia Art and Science, Waseda University, Japan, k.yatabe@asagi.waseda.jp(4)

Department of Intermedia Art and Science, Waseda University, Japan, yoikawa@waseda.jp

Abstract
Visualizing sound fields is important to understand them intuitively. Recently, sound field visualization using aparallel phase-shifting interferometer (PPSI) has been proposed. This optical method can observe a sound fieldinstantaneously and quantitatively without placing any object inside the field. Thus, sound fields difficult forthe ordinary instruments to measure can be investigated by PPSI such as the fields inside a small cavity or air
flow. After a measurement, the observed data must be analyzed to obtain meaningful information in terms ofacoustics. However, such analysis has not been studied much as PPSI itself is a newly developed method. In
this paper, we estimate sound intensity from the data captured by PPSI. The number of the observation pointsof our PPSI system is up to 262,144, where the interval between the adjacent points is 0.22 mm. Therefore,sound intensity can be estimated densely at quite a lot of points. The accuracy of estimated sound intensity isinvestigated through numerical experiments, and then real data observed by PPSI are analyzed for visualizing
the sound intensity.
Keywords: Sound visualization, optical measurement, least squares method, Polynomial approximation

1

INTRODUCTION

Measurement and visualization of sound fields are important because they lead to intuitive understanding. So
far, various measurement and visualization methods have been proposed. For example, a visualization system of
3D acoustic intensity map using mixed reality technology has been proposed [1]. As another example, in recent
years, optical methods [2, 3] for sound fields visualization using parallel phase-shifting interferometer (PPSI)
have been proposed [4–7]. Since changes in air density due to sound change the refractive index and affectlight, it is possible to obtain sound information by using an optical system. By using an optical method, it
is possible to measure from a distance without causing sound diffraction or reflection due to the installationof a measuring instrument. In particular, in the case of a small area near a sound source or a sound field,
due to the relatively large effect of installing a measuring instrument, a non-contact optical method is effective.
Another advantage of the optical method is that calibration is not necessary for quantifying sound pressure
if the wavelength of light and atmospheric pressure are known. Furthermore, PPSI measures multiple points
simultaneously and instantly with a high-speed camera. PPSI can measure the phase difference of light changedby the sound, and the information can be obtained as a movie. The number of the observation points of our
PPSI system is up to 262,144, where the interval between the adjacent points is 0.22 mm, and we can observe
the sound field in a small range densely. Therefore, it is possible to observe time-varying phenomena suchas the sound field of an instrument sound. As an example of visualization of time-varying phenomena, PPSI
was used to simultaneously visualize the flow of air and sound [8, 9], which leads to the understanding of thegeneration mechanism of aerodynamic noise. Therefore, visualization of the sound field using PPSI is useful invarious fields as a technology for understanding acoustic phenomena.
∗ Current affiliation is NTT Comunication Science Laboratries, Nippon Telegraph and Telephone Corp.

2765

In order to apply the values measured by PPSI to other studies, it is necessary to make an acoustic interpretation. In addition to visual interpretation of the observed values as they are, other interpretations can be madeby performing signal processing. However, because research by PPSI has recently started, its application using
signal processing has not been widely conducted compared with other acoustic research, and it is a field underdevelopment [10–12]. Therefore, we aimed to acquire information that matches the characteristics of PPSI, thevisualization of sound fields.
Sound intensity, also known as acoustic intensity, is important information of sound that helps the understanding of the sound field. The sound intensity represents the power of the sound and its energy flow as avector quantity and has information about the direction of the sound. Therefore, even in areas where interpretation is difficult using sound pressure information alone, sound wave propagation can be further understood, and
it is useful in various fields such as sound source localization, sound source separation, architectural sound, andso on. In fact, research has been conducted to make an intuitive understanding of the sound field by visualizing
sound intensity, not sound pressure [1]. Also, because the data of PPSI contains abundant spatial information,
it is reasonable to consider the sound intensity that can be estimated using the spatial information.
In this paper, we propose a method to estimate acoustic intensity using polynomial approximation fromacoustic data captured by PPSI. The conventional method uses the first-order difference approximation. On theother hand, we estimate the intensity using polynomial approximation with second or higher order and using
relatively many observation points. Then, the estimation performance can be expected to be improved by theincrease of order using polynomial approximation and the increase of spatial information with multipoint observation. We perform numerical simulations to apply the proposed method to an ideal point source sound fieldwith noise and compare the performance with that of the conventional methods using frequency and number ofobservation points as variables.

2

SOUND VISUALIZATION BY PPSI

2.1 Principle of sound visualization by interferometer
Acousto-optic interaction is characterized by the sound changing the refractive index of the medium. The relationship between the phase of light and refractive index can be written by a geometrical optics approximation:∫
ϕ(r, t) = kn(l, t)dl,
(1)
L(r)

where ϕ is the phase of light, r is the position vector, t is the time, k is the wave number of the light, L isthe optical path and n is the refractive index. The relationship between the refractive index and density of gas
can be written by the Gladstone–Dale relation [13]:n−1= const.,ρ(2)where ρ is the density. Assuming that the density change of the medium due to sound is an adiabatic change,
the relationship between density and sound pressure follows
p0 + p
=
p0

(

ρ
ρ0

) γ1
,

(3)

where p0 is the pressure in static state, p is the sound pressure, ρ0 is the density in static state, and γ is the
specific heat ratio. From Eqs. (2) and (3), the relationship between the refractive index and the sound pressure
follows
(
)1
p γ
+ 1.
(4)
n = (n0 − 1) 1 +
p0

2766

When p is much smaller than p0 , Eq. (4) can be linearized asn ≃ n0 +n0 − 1p.
γp0(5)From Eqs. (1) and (5), except for steady-state terms that do not depend on sound pressure, the relationshipbetween the phase of light and the sound pressure follows∫
n0 − 1ϕ(r, t) = kp(l, t)dl.
(6)
γp0
L(r)
Therefore, the sound pressure can be obtained by observing the phase of light.
2.2 PPSIInterferometry is a method of measuring the phase difference using light interference, and one of them isphase-shifting interferometry that acquires multiple phase-shifted interferograms [14]. In particular, as a spatial
phase-shifting method, PPSI that can acquire multiple images in one shot has been proposed. The schematic
diagram of a PPSI system is shown in Fig. 1. The laser emits monochromatic light, which is split by the
Wollaston prism into two orthogonal linear polarizations. Then, the lights are divided by the beam splitter into
reference light and object light. The reference light is reflected by the reference surface of the optical flat.
The object light passes through the reference surface and through the measurement area. It is reflected by therefractive surface, passes again through the measurement area. The two polarizations are combined and reflectedby the beam splitter. Then, the high-speed camera observes light whose polarization state has been changed asinterference fringes.
test section

reflective
surfacebeam splitterreference
surface

：lens

laser

Wollaston prism
spatial filter
quarter wave platehigh-speed polarization camera

Figure 1. The schematic diagram of PPSI system.

3

CONVENTIONAL METHOD

Intensity is defined as the time average of the product of pressure and particle velocity:
I(t) = ⟨p(t)u(t)⟩t ,(7)

where ⟨·⟩t is time average. Intensity is originally a complex value, but the active intensity is generally used assound intensity. Therefore, only the activity intensity is considered.
3.1 p-p method
In the general case, particle velocity is determined using two or more pressures, not directly measured. Such amethod is called p-p method [15]. This method uses a finite-difference approximation of the pressure gradient2767

to determine the particle velocity. Then, the particle velocity in the time domain follows∫
1 t p2 (τ ) − p1 (τ )dτ,u(t) = −
ρ −∞
∆r

(8)

where p1 and p2 are the sound pressure from the two observation points and ∆r is the distance with the twopoints. The average of the pressure between two observation points is used for pressure:p(t) =

p1 (t) + p2 (t)
.
2

(9)

Therefore, the intensity is determined as
p1 (t) + p2 (t)I(t) = −2∫ tp2 (τ ) − p1 (τ )
dτ.
ρ∆r
−∞

(10)

While Eq. (10) is in the time domain, the equation in the frequency domain follows
I(ω) = −1
Im[F1 (ω)F2 (ω)],
ρω∆r

(11)

where ω is the angular frequency, Im[·] is the imaginary part, F1 (f ) and F2 (f ) are the Fourier transformsof the pressures p1 (t) and p2 (t), and over-line indicates complex conjugate. Then, this equation integrated infrequency bands f1 to f2 follows∫ f2
Im[F1 (f )F2 (f )]
1
df.
(12)
I(f1 , f2 ) = −2πρ∆r f1
f
The method of calculating intensity using this equation is called the cross-spectral method and is equivalent to
the intensity calculated using equation Eq. (10) when integrated over the frequency band of the observed soundfield. In general, not only in the case of two observation points, Eq. (11) can be represented asI(ω) = −1
Im[pc ∇p],
ρω

(13)

where pc is the sound pressure at the center of observation points, and ∇p is the sound pressure gradient.
3.2 PAGE method
Recently, Thomas et al. proposed a new intensity estimation method and named the phase and amplitude gradient estimation (PAGE) method [16, 17]. This method is based on the work of Mann et al. [18] that determined
the intensity asI(ω) =

1 2
P ∇ϕ,
ρω c

(14)

where Pc is the amplitude of pc .

4

PROPODED METHOD

The conventional methods, both the p-p method and PAGE method, use difference approximation using a firstorder polynomial in Eqs. (13) and (14) to estimate the gradients. We extend them to higher order polynomial
approximation. The N th-order polynomial model at position vector r = [x, y] follows
h(x, y) =

N N
−n
∑
∑
n=0 m=02768

cm,n xm y n ,(15)where c is a coefficient. Then, polynomial approximation can be obtained by solving the least squares problem
for the observed data g asK1∑
arg min|gk − h(xk , yk )|2 ,(16)
2h(x,y)
k=0

where K is the number of observation points.
This method can extend both Eqs. (13) and (14) by substituting pressure or phase for g in Eq. (16). The
extension of Eq. (13) is written asjI(ω) =
p0 ∇hp (x, y),(17)ρωwhere hp (x, y) is h(x, y) when g = p. Also, the extended Eq. (14) is written as
I(ω) =jp0 ∇hϕ (x, y),
ρω

(18)

where hϕ (x, y) is h(x, y) when g = ϕ.

5

EXPERIMENTS

The performance of the proposed method was evaluated by numerical experiments. Intensity estimation accuracy
is compared for angle and magnitude with different methods and frequencies. We used observation points in thesquare region. The target sound field is the ideal point source sound field with a noise of 10 dB. The frequency
was varied from 2 kHz to 22 kHz at 4 kHz intervals, and the sampling frequency was 44.1 kHz. The incident
angle was 45 degrees to the estimated point, and the distance from the sound source was 31.11 mm (equivalent
to 100 adjacent point intervals). We used two values to evaluate the accuracy of the methods: the error in theangleIˆ · Iθerror = arccos,(19)
ˆ
|I||I|and the error of magnitude|I|%error = 100

ˆ − |I||
||I|
,
|I|

(20)

where the hat is used to represent estimated value.
5.1 Comparison of conventional method and proposed method
Figure 2 shows the experimental results of comparison by the methods with Eq. (13). The cross spectrum
method and the proposed method with 3 × 3 channel result in almost the same for both angle error and estimate
error. On the other hand, the proposed method with 25 × 25 channel results in less error than them, regardlessof frequency, the angle error is less than 3 degree, the estimate error is less than 0.07%
Figure 3 shows the experimental results of comparison by the methods with Eq. (14). The proposed methodwith 25 × 25 channel results in almost less than the other methods. however, the errors at 18 kHz are high.
5.2 Comparison by frequency with proposed method
Figure 4 shows the experimental results of comparison by the frequency with the proposed method. Regardlessof frequency, the performance of the proposed method using Eq. (13) increases as the number of observation
points used for polynomial approximation increases. In the case of 6 kHz and 18 kHz sound field, the estimatederror and the squared error of the proposed method using Eq. (14) are large than that in other cases. Therefore,
it was found that the estimated value might be worse depending on the frequency band.

2769

Figure 2. The intensity estimation error of the methods using the pressure gradient as a function of frequency.

Figure 3. The intensity estimation error of the methods using the phase gradient as a function of frequency.

Figure 4. The intensity estimation error with the proposed method as functions of frequency and the numberof channels. The left figure shows the results with Eq. (13). The center figure shows the results with Eq. (14)and the right figure shows the squared errors at that time.

2770

5.3 Applying proposed method to real date
We apply the proposed method with Eq. (13) to estimate the intensity of the sound field of the flute. The
measuring area is 352 mm × 563 mm and Fig. 5 shows comparing the results of proposed, cross spectrum
and PAGE method. Only the intensity vector in (b) of Fig. 5 shows sound propagation from the lower left
sound source. Therefore, the estimation performance of the proposed method is better than that of conventional
methods. However, some vectors in (b) Fig. 5 point in the opposite direction to the surrounding vector. One of
the causes is that the observation data is too noisy.

0.03

70

Intensity [dB]0(b)Phase [rad]

(a)

-0.03

(c)

(d)

40
Figure 5. Reference image (a) of experiment and intensity of sound field of the flute with proposed (b), crossspectrum (c) and PAGE method (d).

6

CONCLUSIONSIn this paper, we proposed a method using a polynomial approximation to estimate sound intensity from phasedata obtained from PPSI. The proposed method is an extension of the conventional method and can be applied
to both Eqs. (13) and (14). By increasing the number of observation points used for the approximation, we
confirmed that the proposed method is effective for noisy fields. In the case of estimation using phase information of sound, it has been found that the estimated value becomes worse depending on the frequency band. In
future research, in addition to the relationship with the frequency band, it is necessary to consider the differencein performance depending on the shape of the observation range and the number of observation points.

REFERENCES[1] A. Inoue, Y. Ikeda, K. Yatabe and Y. Oikawa, Visualization system for sound field using see-through head
mounted display, Acoust. Sci. and Tech., vol. 40 (1), pp. 1–11, 2019.
[2] K. Nakamura, M. Hirayama and S. Ueha, Measurements of air-borne ultrasound by detecting the modulation in optical refractive index of air, 2002 IEEE Ultrason. Symp., vol. 1, pp. 609–612, 2002.
[3] Y. Oikawa, M. Goto, Y. Ikeda, T. Takizawa and Y. Yamasaki,

Sound field measurements based on re-2771

construction from laser projections, IEEE Int. Conf. Acoust. Speech Signal Process., vol. 4, pp. 661–664,
2005.
[4] K. Ishikawa, K. Yatabe, N. Chitanont, Y. Ikeda, Y. Oikawa, T. Onuma, H. Niwaand and M. Yoshii, Highspeed imaging of sound using parallel phase-shifting interferometry, Opt. Express, vol. 24, (12), pp. 12922–
12932, 2016.
[5] Y. Oikawa, K. Yatabe, K. Ishikawa and Y. Ikeda, Optical sound field measurement and imaging using
laser and high-speed camera, 45th Int. Congr. Expo. Noise Control Eng., pp. 258–266, 2016.
[6] K. Ishikawa, K. Yatabe, Y. Ikeda, Y. Oikawa, T. Onuma, H. Niwa and M. Yoshii, Optical sensing of soundfields: non-contact, quantitative, and single-shot imaging of sound using high-speed polarization camera,
5th Jt. Meet. of Acoust. Soc. Am. and Acoust. Soc. Jpn., Nov. 2016.
[7] K. Ishikawa, K. Yatabe, Y. Ikeda, Y. Oikawa, T. Onuma, H. Niwa and M. Yoshii, Interferometric imagingof acoustical phenomena using high-speed polarization camera and 4-step parallel phase-shifting technique,Proc. SPIE, vol. 10328, 103280I, Nov. 2017.
[8] K. Ishikawa, R. Tanigawa, K. Yatabe, Y. Oikawa, T. Onuma and H. Niwa, Simultaneous imaging of flow
and sound using high-speed parallel phase-shifting interferometry, Optics Letters, vol. 43, (5), pp. 991–994,
Mar. 2018.
[9] R. Tanigawa, K. Ishikawa, K. Yatabe, Y. Oikawa, T. Onuma and H. Niwa, Optical visualization of a fluid
flow via the temperature controlling method, Optics Letters, vol. 43 (14), pp. 3273–3276, July. 2018.
[10] K. Yatabe, R. Tanigawa, K. Ishikawa, Y. Oikawa, Time-directional filtering of wrapped phase for observing transient phenomena with parallel phase-shifting interferometry Optics Express, vol. 26 (11), pp.
13705–13720, May 2018.
[11] Kohei Yatabe, Kenji Ishikawa, Yasuhiro Oikawa, Hyper ellipse fitting in subspace method for phaseshifting interferometry: Practical implementation with automatic pixel selection, Optics Express, vol. 25(23), pp. 29401–29416, Nov. 2017.
[12] Kohei Yatabe, Kenji Ishikawa, Yasuhiro Oikawa, Acousto-optic back-projection: Physical-model-based
sound field reconstruction from optical projections, Journal of Sound and Vibration, vol. 394 (28), pp. 171–
184, Apr. 2017.
[13] J. H. Gladstone and T. P. Dale, Researches on the refraction, dispersion, and sensitiveness of liquids,
Philos. Trans. R. Soc. Lond., vol. 153, pp. 317–343, 1863.
[14] T. Yatagai, B. J. Jackin, A. Ono, K. Kiyohara, M. Noguchi, M. Yoshii, M. Kiyohara, H. Niwa, K. Ikuo
and T. Onuma, Instantaneous phase-shifting Fizeau interferom-etry with high-speed pixelated phase-mask
camera, Proc. SPIE 9660, SPECKLE 2015: IV International Conference on Speckle Metrology, 966018, Aug.
2015.
[15] T. D. Rossing, Springer Handbook of Acoustics, Springer, pp. 1093-1103, 2014.
[16] D. C. Thomas, B. Y. Christensen and K. L. Gee, Phase and amplitude gradient method for the estimationof acoustic vector quantities, J. Acoust. Soc. Am., vol. 137, pp. 3366–3376 2015.
[17] K. L. Gee, T. B. Neilsen, S. D. Sommerfeldt, M. Akamine and K. Okamoto, Experimental validationof acoustic intensity bandwidth extension by phase unwrapping, J. Acoust. Soc. Am., vol. 141, pp. EL357–
EL362, 2017.
[18] J. A, Mann III, J. Tichy, and A. J. Romano, Instantaneous and timeaveraged energy transfer in acoustic
fields, J. Acoust. Soc. Am., vol. 82, pp. 17—30, 1987.

2772