**3. Geometrical optical method for modeling aero-optical transmission**

#### **3.1 Ray tracing model**

Geometrical optics is used in this chapter because wavelengths are considered to be negligible. According to the principle of geometrical optics, light exists in the form of straight lines in a uniform medium. When light passes through two homogeneous media with different refractive indices, its behavior can be determined by the laws of refraction and reflection.

According to the Gladstone-Dale relationship above, CFD grids can be converted to indexed grids. The beam axis is oriented parallel to the negative Z direction. As geometric optics show, the refracted/reflected rays are in the same plane as the incident rays. Therefore, the transport of light in a 3D flow field can be seen as consisting of the transport of multiple layers in a 2D cross-section of the flow field. Light incident on CFD grid points now starts at the top of the computer field. **Figure 5** shows how light travels in a plane.

**Figure 5.** *Geometry of the optical transmission in two dimensions.*

The index of refraction at grid point 1–1 is denoted by *n*<sup>11</sup> (see **Figure 3**). Likewise, the index of refraction at point i-j is denoted by *nij*. Points 1–1 to 1–2 are chosen as *n*<sup>11</sup> for the refractive index of the region above the interface, and the index in the next grid of the threshold is selected as *n*21. *d* is the size of the square grid taken from a hexahedral grid plane. Assuming that the initial angle of incidence is *θ*0, the angle of refraction of the ray as it passes through the interface is*θ*1. Thus, when a ray passes through the grid in the negative Z direction, the angle of refraction at point k is written as*θk*. The ray coordinate is (*xk, zk*) and its offset is denoted by Δ*xk*. The offsets are expressed byΔ*x*<sup>1</sup> ¼ j j *X*1*O* ¼ *x*<sup>1</sup> � *O*, Δ*x*<sup>2</sup> ¼ j j *X*2*X*<sup>1</sup> ¼ *x*<sup>2</sup> � *x*<sup>1</sup> … Δ*xk* ¼ j j *XkXk*�<sup>1</sup> ¼ *xk* � *xk*�1. Here, the total real offset is marked byΣΔ*xk*. Through the geometrical relation, the angle of incidence at point 2 is equivalent to the angle of refraction at point 1. According to the Snell's law, an equation at grid point 1–1 *n*<sup>11</sup> sin *θ*<sup>0</sup> ¼ *n*21sin*θ*<sup>1</sup> is acquired. At point 2<sup>0</sup> , an equation *n*<sup>31</sup> cos *θ*<sup>2</sup> ¼ *n*<sup>32</sup> sin *θ*<sup>0</sup> <sup>2</sup> is obtained. The point where the transmission is similar to point 2<sup>0</sup> is denoted by *l* 0 . The optical path length (OPL) from point 1 to point 2 is indicated as *OPL*1. Likewise, the length of the path from point k to point k + 1 is denoted by*OPLk*. Total reflection can occur when light is transmitted from an optically denser layer to an optically thinner layer due to the introduction of an interface. But interfaces obviously do not apparently exist in real airflow. Total reflection causes light to rotate. This definitely increases the complexity of the algorithm.

Rayleigh pointed out that the wavefront cannot be changed when the wavefront error between the real and reference wavefronts was less than a quarter wavelength. The refractive indices of two adjacent grids are denoted *ni* and *nj*, respectively. According to the Rayleigh criterion, the wavefront error is negligible if the following equation is satisfied, that is, if the light transmittance is not deformed. That is, the initial light path can be considered unchanged. Rayleigh's criterion is applied to the algorithm. Resort to judging the criterion when total reflection occurs. If yes, the hypothetical interface is considered nonexistent and the light propagates in the same direction. Otherwise, the virtual interface is assumed to exist and full reflection occurs.

$$\left| n\_{\circ} - n\_{i} \right| s < \lambda/4 \tag{4}$$

where *s* is the geometrical path length and 1≤*s* ≤ ffiffi 2 <sup>p</sup> (mm) The results show that the approach is useful and enough for simplifying the calculations. Modeling the propagation through the CFD grids, we can derive the recursive algorithm for tracing the light rays. First of all, the relationship at point 1 is expressed by Eq. (5).

$$\begin{cases} \Delta \mathbf{x}\_1 = d \tan \theta\_1 \\ n\_{21} \sin \theta\_1 = n\_{11} \sin \theta\_0 \end{cases} \tag{5}$$

At point *k*, the ideal relation of the light transmission was gained as

$$\begin{cases} \Delta \mathbf{x}\_k = d \tan \theta\_k \\ n\_{k+1,l} \sin \theta\_k = n\_{k,l} \sin \theta\_{k-1} & (l, k = 1, 2, \dots) \\ OPL\_k = n\_{k+1,l} d/\cos \theta\_k \end{cases} \tag{6}$$

If ΣΔ*xk* >*l* � *d*, where *l* ∈ *N*and *l* is counter, we should modify the offset Δ*Xk* and the *OPLk*. Through the triangular relation, Eq. (7) is obtained by

*Perspective Chapter: Computational Modeling for Predicting the Optical Distortions… DOI: http://dx.doi.org/10.5772/intechopen.106591*

**Figure 6.** *Interpolated mesh cell.*

$$M = \left(\Sigma \Delta \mathbf{x}\_k - l \times d\right) \cot \theta\_k \tag{7}$$

The modified offset and OPL would be acquired by

$$\begin{cases} \Delta \mathbf{x}\_k = M \cot \theta\_k' + (d - M) \tan \theta\_k \\ \text{OPL}\_k = n\_{k+1,l} (d - M) / \cos \theta\_k + n\_{k+1,l+1} M / \sin \theta\_k' \end{cases} \tag{8}$$

If no reflection over the boundary is produced, the relation between *θ<sup>k</sup>* and *θ*<sup>0</sup> *k* should be given by

$$n\_{k+1,l}\cos\theta\_k = n\_{k+1,l+1}\sin\theta'\tag{9}$$

$$
\theta\_{k+1} = \mathbf{90}^\* - \theta\_k' \tag{10}
$$

In the case of total reflection, it is assumed that the transmission direction has no changes if Eq. (4) is satisfied. Then, Eq. (9) should be changed into Eq. (11).

$$
\theta\_k = \theta\_k'\tag{11}
$$

Integrated with Eqs (4)–(11), the recursive algorithm for tracing the ray based on the CFD grids is derived. For the fine resolution, a method to interpolate the discrete flow field data is shown in **Figure 6**.

The index at point 1 is interpolated by *n*<sup>1</sup> ¼ ð Þ *nA* þ *nD =*2.

The index at point 2 is gained by *n*<sup>2</sup> ¼ ð Þ *nA* þ *nB* þ *nC* þ *nD =*4.

The index at point 3 is expressed by *n*<sup>3</sup> ¼ ð Þ *nA* þ *nB* þ *nC* þ *nD* þ *nE* þ *nF* þ *nG* þ *nH =*8. The rest can be inferred by analogy and higher resolution indexed fields are

obtained. **Figure 7b** shows the interpolated index field. This method helps to make the data more contiguous with the initial index data and to adopt the algorithm described above.

### **3.2 Aero-optical analysis**

The aero-optical quantities, measured and calculated, mainly are the wavefronts' distortion, Strehl ratio, and the line of sight error. There is a general assumption that the mean flow fields produce time-averaged blurring and the line of sight error, whereas the turbulence produces jitter and blurring. In the sight of geometrical optics,

**Figure 7.** *(a) The initial index data (b) the interpolated index data.*

the wavefronts' aberrations arise from the OPL changing. OPL is expressed by Eq. (12).

$$OPL = \int\_{ry} ndl \tag{12}$$

As the rays penetrate the disturbed air, they pick up the absolute *OPD* as follows11:

$$OPD = \int\_{0}^{L} (n - 1)dl = K\_{GD} \int\_{0}^{L} \rho dl \tag{13}$$

where *L* is the total geometrical path length, and *ρ* is the averaged-density value of the local flow. Additionally, if the central ray is considered as the reference ray, whose OPL is marked *OPLref* , relative OPD will be obtained by Eq. (14) [26].

$$\text{OPD} = \text{OPL} - \text{OPL}\_{\text{ref}} \tag{14}$$

Then, OPD data acquired can be transformed to wavefront phase distortion by using the following formula:

$$
\Delta\phi\left(\overrightarrow{r},t\right) = k\text{OPD} \tag{15}
$$

where *k* is the wave number, and *k* ¼ 2*π=λ*. Therefore, optical path differences directly reflecting the variations of the wavefront phase errors. Another quantity of wavefront aberrations is the root-mean-squared (RMS) optical path difference denoted by *σ*2, that is, it is wavefront variance gained as follows [27].

$$
\sigma^2 = 2K\_{GD}^2 \int\_0^L (\rho')^2 l' dl\tag{16}
$$

where *ρ*<sup>0</sup> is the fluctuation density and *l* 0 is the turbulence length scale calculated by CFD. Wavefront variance is a measure of the dispersion along with the mean optical path length of a wavefront and is important for modeling aero-optical parameters

*Perspective Chapter: Computational Modeling for Predicting the Optical Distortions… DOI: http://dx.doi.org/10.5772/intechopen.106591*

**Figure 8.** *The LOS deviation angle.*

such as Strehl ratios and OTF for the turbulent flows. Only if the RMS optical path difference is not very large, is the time-averaged Strehl ratio for a given wavefront phase shift approximated as [28].

$$\mathcal{SR} = \exp\left[-\delta\_{\phi}^{2}\right] \tag{17}$$

Here, *δ*<sup>2</sup> *<sup>ϕ</sup>* is denoted the wavefront phase variance and it is obtained by the following relationship [29]

$$
\left.\delta\_{\phi}\right\|^{2} = k^{2}\sigma^{2}.\tag{18}
$$

In **Figure 6**, the line of sight error is described in terms of optical beam path reversibility. The LOS errors result in the coordinate position excursion of the image formed by the light propagation through the supersonic flow fields. Generally, *q* is of the small scale. So, the arc length *s* is approximately equal to *q*. If *R* is selected as the radius, the LOS deviation angle denoted by *τ* will be expressed as

$$
\pi = q/\mathbb{R}.\tag{19}
$$

Here, the unit of *τ*is rad. *q* is acquired by Eq. (20)

$$q = h \tan \theta - \Sigma \Delta X\_i. \tag{20}$$

Here, the maximum displacement angle *τ* max representing the maximum displacement of the image position. In **Figure 8**, when R = |OB|, the value of *τ* is smaller than the real angle, and when R = |OA|, the value of *τ* is large. Using each CFD grid as a CCD imaging unit, the detection window is abstracted into an image plane consisting of 64 � 64 pixels. The maximum displacement angle directly reflects the maximum displacement of the intensity distribution, that is, the maximum displacement of the image.

### **3.3 Simulation results**

The supersonic flow field CFD data were calculated using the large-eddy simulation (LES) method. All angles of attack in a flow field simulation are the same. The

**Figure 9.** *Solution of density fluctuation along flow direction. Solution at H = 35, and Ma = 7.*

distribution of density fluctuations in the flow field at a height of 35 km and Mach number 7 are shown in **Figures 9** and **10**. In **Figure 9**, density fluctuations near the detection window along the flow direction become bigger and bigger. The distribution of density fluctuations in the positive Y direction is symmetrical. The RMS OPD distribution obtained from Eq. (18) shows the change in phase shift of the wavefront, which in turn shows the characteristics of the flow field in **Figures 11** and **12**. The RMS

**Figure 11.** *The RMS wavefront errors in the flow direction (H = 35 km, Ma = 7).*

*Perspective Chapter: Computational Modeling for Predicting the Optical Distortions… DOI: http://dx.doi.org/10.5772/intechopen.106591*

**Figure 12.**

*The RMS wavefront errors in the Y positive direction (H = 35 km, Ma = 7).*

optical path difference increases along the flow direction, whereas, in general, in the Y direction, the RMS OPD is the greatest at the center of the window and decreases from the center to the other.

Strehl ratios are shown in **Figures 13** and **14**. On the whole, light intensity is weakened along the flow direction, while in Y positive direction, the light intensity at

**Figure 13.** *Strehl ratio in the flow direction (H = 35 km and Ma = 7).*

**Figure 14.** *Strehl ratio in the Y positive direction (H = 35 km and Ma = 7).*

the center of the window is reduced to minimum and is decreased from here to other two sides. Apparently, all the calculation results are qualitatively correct.

The absolute differences of the optical path in the positive Y direction are shown in **Figures 15**–**17**. It can be seen in **Figure 15** that the optical path difference increases as the angle of incidence increases. Therefore, it can be concluded that the sensor within the detection window needs better incident light to reduce wavefront distortion. It can be seen in **Figure 16** that the optical path difference value decreases as the height increases. It can be seen from the standard atmospheric table that the atmospheric density in the range from 0 to 80 km decreases with increasing altitude. Of course, as the free air flow density becomes thinner, the density near the detection window inevitably becomes thinner at higher altitudes. Although the supersonic flux field near the window is compressible, the change in density is small. Therefore, the aberration of the accumulated optical wavefront is still reduced. **Figure 17** shows the optical path differences at different Mach numbers, at the same height and at the same angle of incidence. In aerodynamics, the effect of the compression ratio becomes greater as the velocity of the flow field increases. As the Mach number increases, the density necessarily increases. Therefore, the optical path difference becomes large.

**Figure 15.** *Comparison of OPD at different angles of incidence (H = 40 km and Ma = 7).*

**Figure 16.** *Comparison of OPD at different altitudes (Ma = 7, and* θ*<sup>0</sup> = 5o ).*

*Perspective Chapter: Computational Modeling for Predicting the Optical Distortions… DOI: http://dx.doi.org/10.5772/intechopen.106591*

#### **Figure 17.**

*Comparison of OPD at different Mach Numbers (H = 30 km, and* θ*<sup>0</sup> = 5o ).*


#### **Table 1.** *The maximum deviation angles.*

The maximum deviation angles calculated at different flow fields are given in **Table 1**. It is magnified with the increasing velocity of the flow and reduced with the altitude being higher.
