# Characterization of rough fractal surfaces from backscattered radar data.

I. INTRODUCTION

The scattering of electromagnetic (EM) waves from rough surfaces has been for decades a very interesting subject for scientific investigation. In many cases the main purpose of this research is the characterization of rough surfaces from scattered EM wave data for remote sensing applications, in the microwave or optical regime [1]--[16]. These surfaces can be modelled mathematically with deterministic or random functions [1]-[3]. However, introducing the fractal geometry, these surfaces can be described in a more detailed way in multiscale [1], [3], [8], [17].

In this paper the scattering of EM waves from rough surfaces using the Kirchhoff approximation is examined [1], [2]. In particular, in Section II the mathematical fundamentals for scattering of EM waves from fractal surfaces are summarized [1]-[3].

In Section III our simulation results for the characterization of the rough fractal surfaces from backscattered EM wave data are presented. Finally, conclusions and future related research of ours are described in Section IV.

II. PROBLEM GEOMETRY AND MATHEMATICAL FORMULATION

The geometry of the problem is shown in Fig. 1. More specifically, an incident EM plane wave illuminates a one--dimensionally rough fractal surface extending from x = - L to x = L, as shown in Fig. 1 (note that a generalization to a two - dimensionally rough fractal surface can be easily accomplished [1], [2]). The angle of incidence of the EM wave is [[theta].sub.i] with respect to the vertical z axis, as shown in Fig. 1, where the incident and scattered wave vectors are denoted by [k.sub.i] and [k.sub.s] respectively [1].

Following [1], and in order to describe the surface roughness, a one-dimensional fractal function is used [1], [3], [4]. This fractal function is described by the following equation

[f.sub.r](x) = [sigma]C [N-1.summation over (n=0)] [(D-1).sup.n] sin([K.sub.0][b.sup.n]x + [[phi].sub.n]), (1)

where D (1 < D < 2) is the fractal dimension of the one dimensional curve generating the fractal surface [1], [K.sub.0] = 2[pi]/[[LAMBDA].sub.0] is the fundamental spatial wavenumber of the fractal surface, [[LAMBDA].sub.0] is the corresponding fundamental spatial wavelength, b (where b > 1) is the spatial frequency scaling parameter, [[phi].sub.n] are arbitrary phases and N is the number of tones describing the surface. The amplitude control factor C is given by

C = [{2d(2-D)/[1-[(D-1).sup.2N]]}.sup.1/2] (2)

so that the above surface function (1), has standard deviation (rms height) equal to a (see [1]). Moreover, it can be easily shown that the above surface function (1), exhibits the property of self-similarity, as described by the following equation [1]

[f.sub.r] (x) [approximately equal o] 1/D-1 [f.sub.r] (bx). (3)

In addition, note that since the spatial frequencies [K.sub.0][b.sup.n] are not integer multiples of each other, the above surface function (1) is not a periodic function, but an 'almost periodic function' [1], [3]. Furthermore, it can be easily realized from (1) above, that when the surface fractal dimension D increases from value 1 to value 2 (i.e. 1 < D < 2), the surface roughness also increases (see [1], [3] for more details).

In order to calculate the scattered field from a rough fractal surface, with problem geometry as described in Fig. 1, the Kirchhoff approximation is used in this paper, for which it is assumed that the wavelength of the incident EM wave is small compared to the local radius of curvature of the surface roughness [1]-[3]. Furthermore, for the plane EM wave incidence of Fig. 1, in [1] it is shown that the scattered electric field is given by the following equation

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (4)

where:

p = (1 - R)sin [[theta].sub.i] + (1 + R)sin [[theta].sub.s], (5)

q = (1 + R) cos [[theta].sub.s] - (1 - R) cos [[theta].sub.i], (6)

[[upsilon].sub.x] = k(sin [[theta].sub.i] - sin [[theta].sub.s]), (7)

[[upsilon].sub.z] = -k (cos [[theta].sub.i] + cos [[theta].sub.s]). (8)

In the above (4)-(8) [R.sub.0] is the distance from the observation point (monostatic radar) to the origin, coinciding with the 'source surface point', k is the wavenumber of the incident EM wave (k = 2[pi]f/c, where f is the frequency of the incident EM wave), R is the Fresnel reflection coefficient of the tangential plane at the point of interest, [[theta].sub.s] is the direction of the observer and [f'.sub.r] is the derivative with respect to its argument (x). In this paper, for simplicity we assume a perfectly conducting rough surface, in which case the Fresnel reflection coefficient is given by ([R.sup.+] =1, [R.sup.-]=1), where the superscript + indicates the parallel (vertical) polarization and the superscript - denotes the perpendicular (horizontal) polarization, respectively [1], [2].

In the case of a smooth and perfectly conducting surface, the scattered field for horizontal polarization can be found in the direction of specular reflection, namely for [[theta].sub.i] = [[theta].sub.s] [1], [2]

[E.sub.sc0] = i2k[L.sup.2] exp(ik[R.sub.0])cos [[theta].sub.i]/[pi][R.sub.0]. (9)

By normalizing the value of the scattered field [E.sub.sc] of (4) by the value provided by (9), the scattering coefficient [gamma] is calculated by [1]

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (10)

The first term in the parenthesis provides the most significant contribution to the scattering process, while the second term represents an edge effect, which can be considered negligible when L >> [lambda], as assumed in this paper.

Finally, the scattering coefficient can be expressed in a closed-form solution as follows [1]

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (11)

where [J.sub.m] represent Bessel functions of first kind and order m, and function sinc(x) = sin(x)/x. Furthermore, the superscripts + and - denote the parallel (vertical) and perpendicular (horizontal) polarization respectively, as mentioned above. It must be noted here that, in our simulations, (10) is used, in MATLAB environment, thus exhibiting some advantages over the numerical implementation of (11). However, (11), in addition to (10) offers a significant physical insight to the EM scattering problem considered here.

III. SIMULATION RESULTS

In this paper we concentrate on the backscattering of EM waves from rough fractal surfaces (e.g. monostatic SAR radar [9], [10]), i.e. [[theta].sub.s] = -[[theta].sub.i] at Fig. 1 and (5)-(8), and we plot the backscattering coefficient [absolute value of [gamma](k)] (in magnitude). The surface is simulated as a zero-mean, band-limited fractal function, as in (1), and its roughness is controlled by the fractal dimension D [1], [3]. For example, for D = 1.05 the surface is almost sinusoidal [1] and the roughness is gentle, while as parameter D increases, the surface roughness also increases (1 < D < 2), see [1], [3].

In our simulation experiments, the backscattering coefficient [absolute value of [gamma](k)] was calculated from (10) for a number of frequencies [f.sub.m] = [f.sub.o] + (m - 1)[DELTA]f, where m = 1, 2, ..., M and M is the number of frequencies, fo is the carrier frequency, [DELTA]f = B/M is the frequency step and B is the bandwidth of the radar, i.e. 'stepped - frequency' transmitted radar waveform [9], [10]. Once again, in all our simulations it is set [[theta].sub.s] = -[[theta].sub.i], i.e. backscattered EM wave. Furthermore, in Fig. 2-Fig. 5, below, the plots of [absolute value of [gamma](k)] for angle of incidence [[theta].sub.i] = 30[degrees] are shown, while the values of the other parameters are B = 1 GHz (radar bandwidth), [f.sub.0] = 10 GHz (initial radar carrier frequency) and M = 200 (i.e. 200 frequency steps in radar emitted stepped-frequency waveform).

As far as the simulated fractal surface is concerned, the frequency scaling parameter was set equal to b = 1.8 while the number of tones was set equal to N = 6 [1]. Moreover, the rms height of the surface was set equal to [sigma] = 0.05[lambda], [[LAMBDA].sub.0] = 10[lambda] = 0.3 m and the illuminated length of the rough surface along x-direction ('patch size') was chosen to be 2L = 802 (Fig. 1) in all calculations (so as 2L >> [[LAMBDA].sub.0] and k[sigma] < 1), where [lambda] = c/[f.sub.o] [1], [18]. Note that the length 2L, in practical terms, is proportional to the SAR radar vertical beamwidth, for a SAR radar platform moving in the y-direction of Fig. 1 (not shown in this figure), as well as proportional to the distance [R.sub.0] between the radar and the surface scattering center (e.g. proportional to the height of the SAR radar platform).

Furthermore, at the top left corner of each figure a sample plot of the roughness fractal function [f.sub.r] (x) (1) is also shown.

The roughness of the simulated fractal surface (the fractal dimensionD) is increasing per image, e.g. D = 1.05 (Fig. 2), D = 1.30 (Fig. 3), D = 1.55 (Fig. 4), D = 1.80 (Fig. 5).

By observing Fig. 2-Fig. 5, the following conclusion is made: as the value of the parameter D increases, i.e. as the roughness of the fractal surface increases, the emerging slope between the main lobe and the side lobes also increases.

Therefore, it becomes clear in our simulations that the roughness of the fractal surface can be characterized by the mean slope between the main lobe of function [absolute value of [gamma](k)] and the two sidelobes, adjacent to the main lobe (Fig. 2- Fig. 5). Indeed, as D[right arrow]1 (1 < D < 2) we expect the above slope to follow the slope of a 'sinc' function [sinc(x)=sin(x)/x], see (11), while as D increases (1 < D < 2) we expect, by physical intuition, larger slopes of the function [absolute value of [gamma](k)], due to increased rough surface 'randomness' (at least in the case that the wavelength of the incident wave is 'comparable' to the rough surface 'roughness', i.e. in the 'resonance regime').

In Table I, below, the relation between the fractal dimension D and the slope calculated from each graph is shown. The slope is equal to [absolute value of [DELTA][gamma]]/[absolute value of [DELTA]k], where [DELTA][gamma] represents the amplitude difference between the peak of the main lobe and the peak of the first side lobe, while [DELTA]k represents the difference of the wavenumbers where these peaks occur.

If the radar bandwidth decreases, then the information provided by the backscattered signal-wavenumber plots, of the type provided above, is not always enough in order to draw safe conclusions regarding the roughness (fractal dimension) of the surface. In other words, the bandwidth for our proposed method of surface characterization from backscattered radar data must be sufficiently large (at least 5 % of the carrier frequency [f.sub.o]), in order that the information contained in the plots of Fig. 2-Fig. 5 to be observable and measurable. However, the possible problem of small available bandwidth in actual radar measurements can be counter - balanced by increased 'patch size of observation' 2L [1] [see also the argument of the 'sinc' function in (11). Please note that the 'patch size' 2L is proportional to the SAR radar vertical beam width, as well as to the distance [R.sub.0] between the radar and the surface scattering center. In the following simulation (Fig. 6), the backscattering coefficient [absolute value of [gamma](k)] was calculated for the same fractal dimension D = 1.55 and for two different patch sizes, 2L = 40[lambda] and 2L = 120[lambda]. The rest of the parameters used in this simulation are the same as in Fig. 2-Fig. 5. The corresponding results are shown in Fig. 6. As one can observe in (11) [argument of the 'sinc' function], the number of the sidelobes increases with increasing 'patch size' 2L. Then, by increasing 'patch size' 2L, the steeper become the sidelobe slopes, and the more robust becomes the sidelobe slope calculation.

In order to study further the relation between the surface fractal dimension D and the slopes of the scattering coefficient [absolute value of [gamma](k)], described above, some additional simulations have been performed, as follows. The scattering coefficient [absolute value of [gamma](k)] was calculated sequentially for different values of fractal dimension D, namely here for D = 1.05, 1.10, 1.15, ..., 1.85, 1.90 (i.e. here for 18 subsequent values of parameter D), while the rest of the parameters used in these simulations remained the same as in Fig. 2-Fig. 5. The left and right slope calculations of [absolute value of [gamma](k)], corresponding to each value of D, were averaged, thus creating one average slope calculation of [absolute value of [gamma](k)] for each value of D. Furthermore, in order also to demonstrate the robustness of our proposed method, we inserted in (1) a uniformly distributed random phase variable [[phi].sub.n] in the interval [0, 2[pi]], for every new surface simulation run, namely for every different D value. The calculation results, after 10 simulations for each D value, are presented in Fig. 7, below.

Furthermore, by 'inverting' the data (slope calculations and fitting curve) provided in Fig. 7, the plots of Fig. 8 are provided, as shown below, where in this case the surface fractal dimension D is plotted as a function of the 'slope calculations' of the scattering coefficient [absolute value of [gamma](k)].

In Fig. 8 an analytical expression between the 'slope calculation' and the surface fractal dimension D has been calculated numerically by 'curve fitting'. Namely, it was found here that an excellent curve fitting exists if the fractal dimension D follows the relation D = a * [slope.sup.b] + c . In Fig. 8 we plot the fit curve by using the following coefficients: a = 2.29, b = 0.25, c = 0.913, and as a measure of fit to our calculations we use 'R--square criterion for curve fitting', which, for the above simulations, yielded the value [R.sup.2] = 0.9934 ([R.sup.2] = 1 means 'perfect curve fitting'). We plot also the prediction bounds for the fitted curve. The prediction is based on the existing fit to our simulation calculations by using a 90 % 'probability of occurrence'. Based on the fitted curve model and the prediction bounds for the fitted curve, Tables II and III are presented.

Table II presents the value D that was used for the simulation, the slope that was calculated from this simulation and the Dcalc value, which was calculated using our model. For D = 1.05 to D = 1.15, namely for almost smooth surfaces, the slope is almost zero and for these cases our model could not establish a clear Dcalc value (since the surface is very smooth, in this case). However, for rough surfaces with fractal dimension D > 1.25, our proposed model proved to predict the fractal dimension D of the rough surface with excellent accuracy (see Table II).

The accurate results of Table II gave us the motivation to stress the robustness of our method, by adding random phases [[phi].sub.n] to the rough surface modeling function (1), for each surface simulation. The variability of slope calculations for each value of D in Fig. 7, Fig. 8 depicts this added surface randomness. Table III presents the prediction intervals for the [D.sub.calc] estimation. This interval indicates that for any new observation from a fractal surface characterized by a fractal dimension D, there exists 90 % probability that the value of the fractal dimension of the surface lies within the prediction bounds of this Table (Table III). Summarizing on this, this Table demonstrates the fact that our proposed method is accurate and robust enough to characterize a rough fractal surface from backscattered radar data, and also provides a bound estimation for the fractal dimension D of this surface.

Note that if the angle of incidence is [[theta].sub.i] [approximately equal to] [pi]/2 (i.e. EM wave incidence almost parallel to the rough surface, see Fig. 1), and since in this paper we are interested only for the backscatter case, i.e. [[theta].sub.s]= -[[theta].sub.i], from (8) it follows that in this case [[upsilon].sub.z] [approximately equal to] 0, and from (10) we see that our method is not applicable in this special case (i.e. radar not airborne or spaceborne, any more [9], [10]), since in this special case the scattering coefficient [absolute value of [gamma](k)] cannot be computed. Summarizing, our proposed method for characterizing a rough fractal surface provides reliable results for appropriate values of radar bandwidth, surface 'patch size' and angle of incidence.

IV. CONCLUSIONS

In this paper, a novel method is presented for the characterization of rough fractal surfaces from backscattered SAR radar data of sufficient bandwidth [9], [10]. As resulted from the plots of the backscattered signal magnitude as a function of the wavenumber (frequency) of the incident EM wave, as the roughness of the fractal surface increases (i.e. the fractal dimension D increases), then the observed slope between the main lobe and the side lobes also increases. Moreover, the fractal dimension of the surface can be estimated by the average slope of backscattering coefficient Y(k)|. Furthermore, the value of the available radar bandwidth is crucial and must be sufficiently large, for correct rough surface characterization. In addition, possible lack of enough radar bandwidth for the above purposes can be counter - balanced by increased surface 'patch size' (2L), resulting from increased radar vertical beamwidth or increased radar platform height, as explained in detail above.

V. FUTURE RESEARCH

Regarding future related research, this may include: (i) a three-dimensional (3D) rough fractal surface and scattering modelling which represents a relatively easy extension to the research presented here [1], [2], and it is obviously closer to 'real world' rough surface geometry [1], (ii) sea state characterization by using measured SAR radar data, e.g. for a 'stepped-frequency' (SF) SAR radar waveform, as that used in this paper [9], [10].

http://dx.doi.org/10.5755/j01.eie.22.6.17226

Manuscript received 4 February, 2016; accepted 29 June, 2016.

ACKNOWLEDGMENT

The authors would like to express their sincere thanks to Prof. N. Ampilova and Prof. I. Soloviev, Faculty of Mathematics and Mechanics, St. Petersburg State University, Russia, for very interesting discussions and very helpful suggestions and remarks of theirs, which stimulated this research at the various stages of its implementation.

REFERENCES

[1] D. L. Jaggard, X. Sun, "Scattering from fractally corrugated surfaces", Journal of the Optical society of America, vol. 7, no 6, pp. 1131-1139, 1990. [Online]. Available: http://dx.doi.org/10.1364/JOSAA.7.001131

[2] P. Beckmann, A. Spizzichino, The Scattering of electromagnetic waves from rough surfaces. Artech House Inc., 1987.

[3] D. L. Jaggard, A. D. Jaggard, P. Frangos, "Fractal electrodynamics: surfaces and superlattices", Frontiers in Electromagnetics, pp. 1-47, 2000.

[4] F. Berizzi, E. Dalle Mese, G. Pinelli, "One dimensional fractal model of the sea surface", in IEE Proc. Radar Sonar Navig., vol. 146, no 1, 1999, pp. 55-64. [Online]. Available: http://dx.doi.org/10.1049/iprsn:19990259

[5] A. K. Sultan-Salem, G. L. Tyler, "Validity of the Kirchhoff approximation for electromagnetic wave scattering from fractal surfaces", IEEE Trans. Geosc. Rem. Sensing, vol. 42, no. 9, pp. 1860-1870, 2004. [Online]. Available: http://dx.doi.org/10.1109/TGRS.2004.832655

[6] M. F. Chen, A. K. Fung, "A numerical study of the regions of validity of the Kirchhoff and small perturbation rough surface scattering models", Radio Science, vol. 23, pp. 163-170, 1988. [Online]. Available: http://dx.doi.org/10.1029/RS023i002p00163

[7] E. Jakeman, "Scattering by fractals", in Fractals in Physics, pp. 5560, 1986. [Online]. Available: https://doi.org/10.1016/b978-0-44486995-1.50011-1

[8] N. Ampilova, I. Soloviev, "On digital image segmentation based on fractal and multifractal methods", in Conf. Proc. (CEMA 2015), Sofia, Bulgaria, 2015, pp. 14-17.

[9] A. Malamou, A. Karakasiliotis, E. Kallitsis, G. Boultadakis, P. Frangos, "Application of a fully automatic autofocusing algorithm for post--processing of synthetic aperture radar images based on image entropy minimization", Elektronika ir Elektrotechnika, vol. 19, no. 6, pp. 95-98, 2013. [Online]. Available: http://dx.doi.org/10.5755/j01.eee.19.6.4573

[10] A. Malamou, C. Pandis, P. Frangos, P. Stefaneas, A. Karakasiliotis, D. Kodokostas, "Application of the modified fractal signature method for terrain classification from synthetic aperture radar images", Elektronika ir Elektrotechnika, vol. 20, no. 6, pp. 118-121, 2014. [Online]. Available: http://dx.doi.org/10.5755/j01.eee.20.6.7281

[11] S. Savaidis, P. Frangos, D. L. Jaggard, K. Hizanidis, "Scattering from fractally corrugated surfaces: an exact approach", Optics Letters, vol. 20, no. 23, pp. 2357-2359, 1995. [Online]. Available: http://dx.doi.org/10.1364/OL.20.002357

[12] S. Savaidis, P. Frangos, D. L. Jaggard, K. Hizanidis, "Scattering from fractally corrugated surfaces using the extended boundary condition method", Journal of the Optical Society of America, vol. 14, no. 2, pp. 475-485, 1997. [Online]. Available: http://dx.doi.org/10.1364/JOSAA.14.000475

[13] D. L. Jaggard, X. Sun, "Fractal surface scattering: a generalized Rayleigh solution", Journal of Applied Physics, vol. 68, no 11, pp. 5456-5462, 1990. [Online]. Available: http://dx.doi.org/10.1063/1.347002

[14] M. F. Chen, A. K. Fung, "A numerical study of the regions of validity of the Kirchhoff and small-perturbation rough surface scattering models", Radio Science, vol. 23, pp. 163-170, 1988. [Online]. Available: http://dx.doi.org/10.1029/RS023i002p00163

[15] E. I. Thorsos, "The validity of the Kirchhoff approximation for rough surface scattering using a Gaussian roughness spectrum", Journal of Acoustical Society of America, vol. 83, no 11, pp. 78-92, 1988. [Online]. Available: http://dx.doi.org/10.1121/L396188

[16] J. M. Soto-Crespo, M. Nieto-Vesperinas, "Electromagnetic scattering from very rough random surfaces and deep reflection gratings", Journal of the Optical Society of America, vol. 6, pp. 367-384, 1989. [Online]. Available: http://dx.doi.org/ 10.1364/JOSAA.6.000367

[17] D. L. Jaggard, Y. Kim, "Diffraction by bandlimited fractal screens", Journal of the Optical Society of America, vol. 4, pp. 1055-1062, 1987. [Online]. Available: http://dx.doi.org/10.1364/JOSAAA001055

[18] J. Fikioris, Introduction to Antenna Theory and Propagation of Electromagnetic Waves. National Technical University of Athens (NTUA): Athens, Greece, 1982, pp. 155-156.

Apostolos Dimitrios Kotopoulis (1), Anna Malamou (1), George Pouraimis (1), Evangelos Kallitsis (1), Panayiotis Vassilios Frangos (1)

(1) School of Electrical and Computing Engineering, National Technical University of Athens, 9 Iroon Polytechniou St., 15773 Zografou, Athens, Greece pfrangos@central.ntua.gr
```
TABLE I. FRACTAL DIMENSION D AND THE RESULTING SLOPE
CALCULATIONS.

Left slope
calculations

D       [DELTA][gamma]   [DELTA][k]   slope

1.05        0.0000          0.00      0.0000
1.30        0.0044          2.51      0.0018
1.55        0.0285          5.03      0.0057
1.80        0.0883          3.98      0.0222

Right slope
calculations

D       [DELTA][gamma]   [DELTA][k]   slope

1.05        0.0000          2.62      0.0000
1.30        0.0014          2.93      0.0005
1.55        0.0202          3.77      0.0054
1.80        0.0819          3.77      0.0217

TABLE II. ESTIMATION OF [D.sub.calc]
([D.sub.calc] = 2.29 * sloped 1/4 + 0.913).

D       Slope    [D.sub.calc]

1.05    0.0000       none
1.15    0.0000       none
1.25    0.0002       1.19
1.35    0.0011       1.33
1.45    0.0031       1.45
1.55    0.0065       1.56
1.65    0.0112       1.66
1.75    0.0183       1.76
1.85    0.0274       1.84

TABLE III. [D.SUB.CALC] PREDICTION INTERVAL USING
PREDICTION BOUNDS.

[D.sub.calc]     D     [D.sub.calc]
lower                     upper

1.06            1.05       1.16
1.07            1.15       1.17
1.21            1.25       1.31
1.27            1.35       1.39
1.40            1.45       1.50
1.50            1.55       1.59
1.61            1.65       1.70
1.71            1.75       1.81
1.80            1.85       1.89
```
COPYRIGHT 2016 Kaunas University of Technology, Faculty of Telecommunications and Electronics
No portion of this article can be reproduced without the express written permission from the copyright holder.