## Abstract

Gas flow from an under-expanding jet nozzle is measured using synchronized background oriented schlieren and interferometry diagnostics. The gas density distribution is obtained from a shift vector field of the background oriented schlieren and compared to the interferometric data. The comparison makes use of a simple calibration routine and open source Python recipes.

## 1. INTRODUCTION

Many applications of miniature gas jets require characterization of the gas
or plasma density distribution. To this end, much research has been
performed in a wide range of experimental studies spanning from small
space thrusters to laser-generated plasma targets. In laser wakefield
accelerators, for example, precise control over the gas and plasma density
and its profile is critical to effective acceleration [1]. The flagship gas density measurement
diagnostic is interferometry [2]. A
typical setup of such a measurement involves a probe laser, beam
splitters, mirrors, and imaging optics. In complex experimental
configurations, when *in situ* gas or plasma
density measurements are required, it might be quite difficult or even
impossible to set up such an interferometer.

In this work, we benchmark a simple and robust diagnostic, background oriented schlieren (BOS), against interferometry. We made the benchmarking process more challenging by performing measurements of a flowing gas using an under-expanded jet nozzle. As a result, the gas density has the pattern associated with the well-known shock diamonds [3]. Shock diamonds have been the objects of previous BOS measurements [4,5]. In this work, we show a simple BOS calibration, provide an analytical model for refractive index calculation, translate this model into a numerical recipe, and compare the results with the interferometry measurements. Description of the principles of the BOS diagnostic can be found in [6].

## 2. EXPERIMENT

The experimental setup is depicted in Fig. 1. A 2 mm diameter, axisymmetric, under-expanded gas nozzle was placed in one arm of a Mach–Zehnder interferometer. The interferometer used a 532 nm, CW laser with a beam expander to cover the area of interest above the nozzle. The nozzle was imaged onto a CCD camera using a 55 mm macro lens. The camera was a Point Gray model FL3-GE-28S4M-C, with ${{1928}} \times {{1448}}$ pixels sensor. The magnification of the system was 0.54 (6.88 µm of real object per one pixel of the camera). The interference fringes were aligned perpendicular to the gas flow direction.

The BOS diagnostic was positioned perpendicular to the interferometer probe beam, as shown in Fig. 1. The back-illuminated random background was printed on a transparency using a commercial laser printer at 1200 dots per inch (DPI) or about 21 µm per dot resolution. The background was generated using a simple python code (see Appendix A). The random pattern of the background was 400 by 400 elements in a 2 by 2 cm area (50 µm per element). The background was illuminated by a 5 W flashlight powered with a DC power supply. The background was imaged by a camera that is identical to the interferometer CCD camera using a 50 mm objective lens with a 10 mm extension tube for close-up imaging. In additional experiments, a commercial telecentric lens was used for imaging (see below). Magnification of the BOS diagnostic varied between 0.4 and 0.5.

The distance between the gas nozzle and the background was varied between 10 and 20 mm. The nozzle was connected to a compressed air reservoir at 55 PSI. The reservoir was large enough to provide constant air flow for 5 s before the pressure dropped down, and the reservoir had to be refilled. In the middle of this 5 s window, both CCD cameras acquired synchronously triggered images with the same exposure time of 100 µs. Immediately after that, the gas was shut off, and both cameras acquired another image. This second set of images was used as the references for the corresponding diagnostics.

## 3. OPTICAL CALIBRATION AND SENSITIVITY OF BOS

The optical setup for the BOS diagnostic is depicted in Fig. 2. The most important parameter of the diagnostic is the distance ${Z_3}$ between the background and the nozzle. This distance defines the sensitivity—a magnitude of the shift of the background pixels caused by the angular deflection of the beam in the gas flow. In the paraxial approximation, the sensitivity is linearly proportional to ${Z_3}$, as shown below.

In general, any medium with a refractive index gradient can be represented as a collection of lens elements (spherical or cylindrical). Spatial variation of the refractive index leads to deflection of the rays in proportion to gradients of the refractive index [7]. The deflection angle $\theta = {h_1}/F$, where $F$ is the focal length of this element, and ${h_1}$ is point of deflection (somewhere in the gas jet). The deflection angle (see Fig. 2) is

The BOS CCD camera must be focused on the background, leaving the nozzle out of focus. For regular (entocentric) imaging optics, objects closer to the lens appear bigger on the camera, and the corresponding optical calibration must be done in the nozzle plane. In general, the magnification of the out-of-focus objects can be calculated, but such calculations require detailed knowledge of the imaging optics. This task might be hard in the case of a compound lens with unknown locations of the principal points. The calibration of the nozzle plane magnification can be done experimentally, for example, by counting the number of pixels in the nozzle diameter. Under optimal conditions of operation (resolution versus sensitivity trade-off [8]), the nozzle image appears blurred on the CCD camera [Fig. 3(a)]. This problem can be resolved by placing a periodic pattern (like transparent ruler) in the nozzle plane. The spatial period of the pattern provides the necessary calibration.

A telecentric lens can provide a simpler, robust, but more expensive solution for the calibration problem. Telecentric lenses were used in BOS measurements [9] and do not require any additional calibration. All objects, in or out of focus, appear with the same magnification on the CCD camera. More expensive commercial telecentric lenses provide an adjustable magnification option. The major disadvantages of the telecentric lenses are a short working distance and limited field of view. The price and size of the lens increase drastically when these parameters are increased. In the present experimental setup, a lower price range lens (Edmund Optics 88385) was sufficient for imaging the gas flow.

## 4. DETERMINATION OF REFRACTIVE INDEX IN AXISYMMETRIC FLOW

In the geometric optics limit, we derive the relationship between the gas jet density and the pixel shift at the CCD plane. The position ${\bf{s}} = x{\hat{\bf x}} + y{\hat{\bf y}} + z{\hat{\bf z}}$ (coordinates shown in Fig. 2) and wavenumber ${\bf{k}}$ of a ray of light are given by [7,10]

Consider the horizontal principle ray (shown as black in Fig. 2) coming from a point, denoted as $0{\hat{\bf x}} + {h_1}{\hat{\bf y}} - {Z_3}{\hat{\bf z}}$, on the background screen. The unperturbed ray has the wavenumber and trajectory given by ${{\bf{k}}_0} = {k_0}{\hat{\bf x}}$ and ${{\bf{s}}_0}(z) = 0{\hat{\bf x}} + y{\hat{\bf y}} + z{\hat{\bf z}}$, where $y = {h_1}$, and we have changed the independent variable from time $t$ to the unperturbed ray’s $z$ position $z = ct/{n_0}$. The transverse deflection relative to ${{\bf{k}}_0}$ is governed by

The integral is along a cord intersecting the gas jet. With a change of variables in integration from $z$ to $r$, then Eq. (10) becomes the Abel transform for ${n_1}$:

The inverse Abel transform can be carried out for both directions, but the $y$ component of the deflection is of particular interest because it is the result of radial gradients in the density. The radial gradient is given by

The angular deflection in the $y$ direction ${\theta _y}$ is not directly measured, but Eq. (1) shows the relationship between an angular deflection $\theta$ and the transverse displacement at the camera $\Delta {h_2}$. Additionally, the coordinates $x$, $y$, and $r$ are in the object plane of the gas density perturbation, not the plane of the camera where the displacement is measured. The measured displacement $\Delta {h_2}$ is a function of camera plane coordinate ${h_2}(x,y) = My$. Therefore, the Abel transform in the plane of the camera is

The gas density approaches the ambient density ${n_0}$ as $r \to \infty$ and allows for integration of the above expression:

The first integration (over ${h_2}$) in Eq. (14) is the inverse Abel transformation from projection coordinate $y$ to radius $r$. The second integration calculates the refractive index from its gradient. As explained in Section 3, pixels of shift in the BOS camera image must be translated into real length units of the gas jet plane. In this case, the magnification $M$ is already accounted for in calibration of $h$ and $r$ values and must be omitted from the numerical calculation of the Eq. (14) integrals presented below.

## 5. NUMERICAL ANALYSIS OF BOS DATA

A vector field representing the pixels’ shift was constructed using the
Python implementation of OpenPIV [11] from the image pair acquired in the experiment. The images
were saved in 16 bits tagged image file format (TIFF) files and imported
into OpenPIV using the *openpiv.tools.imread*
function. The full frame image was cropped to the region of interest (ROI)
shown approximately by the dashed rectangle in Fig. 3(a).

After selecting the ROI, the initial cross correlation was done using the
*openpiv.process.extended_search_area_piv*
function. The size of the interrogation window was between 8 and 12
pixels, and the overlap was 2 pixels smaller. The quiver plot coordinates
were extracted using the *openpiv.process.get_coordinates* function, the outliers were
removed with the *openpiv.validation.sig2noise_val* function, and the missing
vectors were substituted by the *openpiv.filters.replace_outliers* function. The final vector
shift map [Fig. 3(b)] was scaled
from pixels of the camera to millimeters of the gas jet plane using the
*openpiv.scaling.uniform* function. This
procedure is similar to an example in the OpenPIV tutorial [12]. This part of the analysis is the
most time consuming and might take a few minutes on a desktop
computer.

The vertical ($y$ direction in notations of Fig. 2) component of the shift vectors created a schlieren image of the gas flow [Fig. 3(c)]. The gas flow from under-expanded nozzle produced a periodic structure that is characteristic of shock diamonds, and the gradient of the refractive index was encoded in this image. Recall that the gas flow was along the positive $x$ -axis. The Fig. 3(c) image is anti-symmetric with respect to the horizontal axis that passes through the center of the gas flow. The inverse Abel transformation in Eq. (14) represents the vertical shift component divided by the corresponding vertical projection coordinate $\frac{{\Delta {h_2}(x,{h_2})}}{{{h_2}}}$. This ratio is expected to be a symmetric function of ${h_2}$. To ensure this, it is possible to introduce a new anti-symmetric vertical coordinate $y$ in the range $[- {y_{{\max}}}/2,{y_{{\max}}}/2]$ with the line $y = 0$ at the symmetry axis (in $[x,y]$ plane). The final symmetrical function that is suitable for Abel inversion has the appearance of $\frac{{\Delta {h_2}(y) - \Delta {h_2}(- y)}}{{2y}}$. In case the new $y$ coordinate is exactly zero at the axis (odd number of pixels in the vertical dimension of the ROI), the neighbor value is employed.

We used the Python PyAbel library [13] to transform the projection
coordinates of the refractive index gradient into three-dimensional (3D)
cylindrically symmetric coordinates. The symmetry axis in PyAbel must be
vertical, and the horizontally symmetrical image must be transposed prior
to inversion. We used the basex transformation method [14] and, because the integration in
PyAbel uses pixels of the image, applied the experimental scaling factor
(corrected to the overlap parameter used in OpenPIV), to keep the
numerical values in real physical units. Finally, the second integral of
Eq. (14) was calculated
using the standard Python NumPy *cumtrapz*
function applied to the result of the Abel transform multiplied by the
coordinate $y$ introduced above. The result was divided
by factor from Eq. (14) to
get values of the refractive index change. The noise was reduced by a
Gaussian filter, and the final result is shown in Fig. 4.

## 6. RESULTS AND DISCUSSION

The radial dependence of the gas density in the flow can be calculated from the refractive index value using the Gladstone–Dale relation [15]:

where $\eta (r)$ is the gas number density to be determined, ${\eta _0} = 2.5 \times {10^{19}}\;{\text{c}}{{\text{m}}^{- 3}}$ is the number density of unperturbed gas, and ${n_0} = 1 + 2.3 \times {10^{- 4}}$ is the refractive index of unperturbed air. The value of the refractive index change $n(r) - 1$ calculated from Eq. (14) is shown in Fig. 4. The number gas density in the jet is shown in Fig. 5(a).Gas density measurements using interferometry were extracted from the raw interferograms using standard techniques described elsewhere [16]. Inverse Abel transformation of the interferometry data employed the same approach as we used in the BOS calculations. The cross section of the gas density extracted from the interferometry is presented in Fig. 5(b). Lineouts of BOS and interferometry density profiles 0.3 mm from the gas flow axis are plotted in Fig. 5(c). Since the CCD cameras (Fig. 1) were not lined up, we shifted the horizontal position of the curves to match the phase of the density modulations in Figs. 5(c) and 5(d). BOS measurement lineout in Fig. 5(d) is an example of using a telecentric lens in the BOS diagnostic. In both cases, agreement between the diagnostics is remarkable.

To extract BOS data, we used only vertical components of the shift vector field. The horizontal component can also be used for the density calculation, but the simple recipe and analysis provided here requires a different approach and will be discussed in a forthcoming paper. In fact, the horizontal component of the shift is analogous to an interference pattern with fringes oriented along the flow. Since the main interest is to obtain the radial dependents of the gas density, such orientation is not optimal.

## APPENDIX A: BACKGROUND GENERATION CODE

The Python code below generates a ${{50}} \times {{50}}$ pixels random background pattern and saves it to a background.png file. The final image file must be scaled appropriately prior to printing.

`#!/usr/bin/env python`

`# coding: utf-8`

`import numpy as np`

`import matplotlib.pyplot as plt`

`# number of pixels`

`row = 50`

`col = 50`

`# figure size on the screen`

`plt.rcParams[“figure.figsize”] =` [5,5]

`# random pattern generation`

`image = np.random.randint(2, size = (row,
col))`

`# show result and save image to file`

`plt.imshow(image, ’gray’)`

`plt.axis(’off’)`

`plt.savefig(’background.png’, bbox_inches = ’tight’,
pad_inches = 0)`

## Funding

U.S. Naval Research Laboratory.

## Disclosures

The authors declare no conflicts of interest.

## REFERENCES

**1. **S. Bulanov, N. Naumova, F. Pegoraro, and J. Sakai, “Particle injection into the
wave acceleration phase due to nonlinear wake wave
breaking,” Phys. Rev. E **58**, R5257–R5260
(1998). [CrossRef]

**2. **E. P. Goodwin and J. C. Wyant, *Field Guide to Interferometric Optical
Testing* (SPIE,
2006), Vol. FG10.

**3. **M. L. Norman, K. H. A. Winkler, L. Smarr, and M. D. Smith, “Structure and dynamics of
supersonic jets,” Astron. Astrophys. **113**, 285–302
(1982).

**4. **N. Van Hinsberg and T. Rösgen, “Density measurements using
near-field background-oriented Schlieren,”
Exp. Fluids **55**, 1720
(2014). [CrossRef]

**5. **F. Nicolas, D. Donjat, O. Léon, G. Le Besnerais, F. Champagnat, and F. Micheli, “3D reconstruction of a
compressible flow by synchronized multi-camera BOS,”
Exp. Fluids **58**, 46
(2017). [CrossRef]

**6. **M. Raffel, “Background-oriented Schlieren
(BOS) techniques,” Exp. Fluids **56**, 60 (2015). [CrossRef]

**7. **L. D. Landau and E. M. Lifshitz, *The Classical Theory of
Fields*, Vol. 2 of Course of
Theoretical Physics Series
(Butterworth-Heinemann,
1980).

**8. **A. B. Gojani, B. Kamishi, and S. Obayashi, “Measurement sensitivity and
resolution for background oriented Schlieren during image
recording,” J. Visualization **16**, 201–207
(2013). [CrossRef]

**9. **M. Ota, F. Leopold, R. Noda, and K. Maeno, “Improvement in spatial
resolution of background-oriented Schlieren technique by introducing a
telecentric optical system and its application to supersonic
flow,” Exp. Fluids **56**, 48 (2015). [CrossRef]

**10. **T. Khan and A. Thomas, “On derivation of the radiative
transfer equation and its spherical harmonics approximation for
scattering media with spatially varying refractive indices,”
Technical Report
(Clemson University Mathematical
Sciences, 2003),
pp. 1–49.

**11. **https://github.com/OpenPIV/openpiv-python.

**12. **https://openpiv.readthedocs.io/.

**13. **https://zenodo.org/record/3370021.

**14. **D. D. Hickstein, S. T. Gibson, R. Yurchak, D. D. Das, and M. Ryazanov, “A direct comparison of
high-speed methods for the numerical Abel transform,”
Rev. Sci. Instrum. **90**,
065115 (2019). [CrossRef]

**15. **J. H. Gladstone and T. P. Dale, “Researches on the refraction,
dispersion, and sensitiveness of liquids,”
Philos. Trans. R. Soc. London **153**,
317–343 (1863). [CrossRef]

**16. **S. Sugawara, S. Nakao, Y. Miyazato, Y. Ishino, and K. Miki, “Three-dimensional
reconstruction of a microjet with a mach disk by Mach–Zehnder
interferometers,” J. Fluid Mech. **893**, A25 (2020). [CrossRef]