ABSTRACT Title of Document: STUDIES OF THE EFFECTS OF ATMOSPHERIC TURBULENCE ON FREE SPACE OPTICAL COMMUNICATIONS. Heba Yuksel, Doctor of Philosophy, 2005 Directed By: Professor Christopher C. Davis Department of Electrical and Computer Engineering Even after several decades of study, inconsistencies remain in the application of atmospheric turbulence theories to experimental systems, and the demonstration of acceptable correlations with experimental results. This dissertation shows a flexible empirical approach for improving link performance through image analysis of intensity scintillation patterns coupled with frame aperture averaging on a free space optical (FSO) communication link. Aperture averaging is the effect of the receiver size on the power variance seen at the receiver. A receiver must be large enough to collect sufficient power and reduce scintillation effects at a given range, but must also be of practical size. An imaging system for measuring the effects of atmospheric turbulence and obscuration on FSO links will be presented. Weak and intermediate turbulence results will be shown for an 863 meter link at the University of Maryland. Atmospheric turbulence has a significant impact on the quality of a laser beam propagating through the atmosphere over long distances. Turbulence causes intensity scintillation and beam wander from propagation through turbulent eddies of varying sizes and refractive index. This can severely impair the operation of target designation and FSO communications systems. A new geometrical model to assess the effects of turbulence on laser beam propagation in such applications will be presented. The atmosphere along the laser beam propagation path is modeled as a spatial distribution of spherical bubbles with refractive index discontinuity statistically distributed according to various models. For each statistical representation of the atmosphere, the path of rays is analyzed using geometrical optics. These Monte Carlo techniques can assess beam wander, phase shifts and aperture averaging effects at the receiver. An effective C n 2 can be determined by correlating beam wander behavior with the path length. In addition, efficient computational techniques have been developed for various correlation functions that are important in assessing the effects of turbulence. The Monte Carlo simulations are compared with the predictions of wave theory. This is the first report to present weak and intermediate turbulence results using an efficient imaging technique. It is also the first report to geometrically simulate aperture averaging. STUDIES OF THE EFFECTS OF ATMOSPHERIC TURBULENCE ON FREE SPACE OPTICAL COMMUNICATIONS. By Heba Yuksel Dissertation submitted to the Faculty of the Graduate School of the University of Maryland, College Park, in partial fulfillment of the requirements for the degree of Doctor of Philosophy 2005 Advisory Committee: Professor Christopher C. Davis, Chair Assistant Professor Thomas E. Murphy Professor Victor Granatstein Professor Robert Newcomb Professor Robert Gammon ? Copyright by Heba Yuksel 2005 Dedication To my parents for their unconditional love and support AND To my husband and children for giving me the drive to succeed ii Acknowledgements I gratefully thank my advisor, Professor Christopher C. Davis, for his continuous guidance, friendship, and support throughout my graduate studies at the University of Maryland. Dr. Davis has been a true inspiration through his kindness and patience that provides an ideal atmosphere to innovate. I also want to thank the members of the dissertation committee for their time and service. I want to particularly thank Assistant Professor Thomas E. Murphy for his helpful technical discussions that allowed me to rethink and solve some of the technical issues that I faced. I am thankful to all of my colleagues in the Maryland Optics Group at the University of Maryland for the fruitful discussions we have held and assistance we have exchanged throughout my years of graduate study at the University of Maryland. I am particularly thankful to Dr. Igor Smolyaninov, Dr. Quirino Balzano, and Dr. Vildana Hodzic. I am also very much thankful to Shawn Ho, Sugi Trisno, Jaime Llorca, and Clint Edwards for the technical and non-technical discussion we continuously held. I am also thankful to Dr. Stuart Milner for continuously initiating new ideas for me to explore as well as collaborating with me in publications. I want to thank Mr. Jay A. Pyle for building all of the necessary mechanical parts for my experimental setup. Mr. Jay always built the parts in a superior turnaround time that always allowed me to meet my deadlines and I am extremely grateful for that. I also want to thank my undergraduate assistants for all of their efforts in collecting data. I especially want to thank Rouzbeh Azarbayejani for setting up the iii weather station and generating a website to monitor it. I am also very grateful to John Lin for devoting extra time to complete some necessary experimental measurements. I am deeply grateful to Professor Howard Milchberg for introducing me to research at an early age and developing my interest in pursuing graduate studies. I am also extremely thankful to Professor Wesley Lawson who has been my academic advisor since the beginning of my undergraduate studies and for his continuous guidance till today. Finally, I want to thank my family. I can not thank enough my parents, Dr. Mohamed Fathy El-Erian and Magda Elkhawalka and sister, Hala Elerian Flores, for their continuous love, guidance, and support throughout my entire educational path. I truly appreciate my parents? complete devotion in helping me to always progress and succeed and providing all of the necessary support to ensure it happening. I could not have done it without them. I want to thank my husband, Dr. Cuneyt Yuksel, for initiating my drive to pursue my Ph.D. and providing me with all of the necessary emotional support. In addition, I greatly appreciate his outstanding patience for bearing my long years of study and continuous travel from our home in Turkey to complete my research at the University of Maryland. I also thank my children, Mehmet Edip & Yasemin Yuksel, for being a continuous inspiration for me to succeed. iv Table of Contents Dedication..................................................................................................................... ii Acknowledgements......................................................................................................iii Table of Contents.......................................................................................................... v List of Tables .............................................................................................................. vii List of Figures............................................................................................................viii List of Publications ..................................................................................................... xii Chapter 1: Introduction................................................................................................. 1 1.1 Thesis Contributions ..................................................................................... 1 1.2 Turbulence Overview.................................................................................... 2 1.3 Aperture Averaging ...................................................................................... 4 1.4 Geometrical Simulations............................................................................... 6 1.5 Thesis Organization ...................................................................................... 8 Chapter 2: Turbulence Theory Overview ................................................................... 10 2.1 Introduction................................................................................................. 10 2.2 Atmospheric Turbulence............................................................................. 10 2.3 Power Spectrum Models for Refractive index fluctuations........................ 13 2.3.1 Kolmogorov Spectrum........................................................................ 15 2.3.2 Tatarski Spectrum ............................................................................... 16 2.3.3 Von Karman Spectrum ....................................................................... 16 2.3.4 Modified Atmospheric Spectrum........................................................ 16 2.4 Weak Turbulence Theory ........................................................................... 17 2.4.1 The Rytov Approximation .................................................................. 17 2.4.2 Numerical Calculation of Fante?s Correlation Functions in Weak Turbulence .......................................................................................................... 19 2.5 Strong Turbulence Theory .......................................................................... 24 2.5.1 Andrews-Prokhorov Asymptotic Analysis ......................................... 24 2.5.2 Churnside Asymptotic Analysis ......................................................... 25 2.5.3 Numerical Calculation of Fante?s Correlation Functions in Strong Turbulence .......................................................................................................... 26 2.6 Weak to Strong Turbulence Theory............................................................ 29 2.6.1 Scintillation Index Model ................................................................... 33 2.7 Aperture Averaging .................................................................................... 34 2.7.1 Aperture Averaging calculations using Fante?s correlation functions in weak turbulence .................................................................................................. 36 2.7.2 Aperture Averaging calculations using Fante?s correlation function in strong turbulence................................................................................................. 38 2.8 Conclusions................................................................................................. 40 Chapter 3: Experimental Setup and Results................................................................ 41 3.1 Introduction................................................................................................. 41 3.2 Aperture Averaging Experimental Setup and Methodology ...................... 43 3.2.1 Experimental Setup............................................................................. 43 3.2.2 Frame Analysis ................................................................................... 46 v 3.3 Aperture Averaging Results........................................................................ 49 3.3.1 Weak Turbulence Results ................................................................... 52 3.3.2 Intermediate Turbulence Results ........................................................ 54 3.4 Optimizing Aperture Size in Receiver Design............................................ 62 3.5 Shape Independence in Image Frame Analysis .......................................... 62 3.6 Weather Effect on Turbulence Level.......................................................... 65 3.7 Conclusions................................................................................................. 66 Chapter 4: Aperture Averaging Effect on Performance of FSO Links....................... 67 4.1 Introduction................................................................................................. 67 4.2 The Normal Distribution of Log Amplitude Fluctuations .......................... 67 4.3 Bit Error Rates in turbulent FSO links........................................................ 69 4.4 Effective Signal-to-noise Ratio analysis..................................................... 74 4.5 Conclusions................................................................................................. 79 Chapter 5: Geometrical Simulation Models ............................................................... 81 5.1 Introduction................................................................................................. 81 5.2 Theory Overview ........................................................................................ 83 5.2.1 Various Approaches............................................................................ 83 5.2.2 Phase Structure Function .................................................................... 84 5.2.3 Beam Wander...................................................................................... 85 5.3 The Spherical/Bubble Model ...................................................................... 87 5.3.1 Beam Wander Simulation................................................................... 88 5.3.1.1 Simulation Procedure and Geometrical Analysis ........................... 89 5.3.1.2 Beam Wander Simulation Results .................................................. 95 5.3.2 Phase Wander Simulation Model........................................................ 97 5.3.2.1 Phase Wander Simulation Results ................................................ 101 5.3.3 Aperture Averaging Simulation Model ............................................ 102 5.3.3.1 Simulation Procedure and Geometrical Analysis ............................. 103 5.3.3.2 Simulation Results ............................................................................ 108 5.4 Conclusions.............................................................................................. 117 Chapter 6: Conclusions and Future Research ........................................................... 119 REFERENCES ......................................................................................................... 123 vi List of Tables Table 3.1: List of the outputs of the LabVIEW program and their description. ??51 Table 5.1: C n 2 for varying n ? values. A n ? of 1e-7 already gives a relatively weak turbulence level while a n ? of 1e-5 represents relatively strong turbulence. ??..109 vii List of Figures Fig. 2.1: A pictorial description of the process of turbulent decay. As turbulent eddies subdivide, they become smaller and more uniform until all of their energy dissipates [4]. ???????????????????????????.11 Fig. 2.2: Spectrum of the refractive index fluctuation. The energy input range, inertial subrange, and energy dissipation ranges are indicated [11]. ??????..15 Fig. 2.3: Von Karman spectrum versus the wave number, ?. ?????????.22 Fig. 2.4: Correlation Function ),( ? ? LB versus the separation between the points ? in weak turbulence using our link parameters. ???????????????...23 Fig. 2.6: Fante?s g(W) function plotted using the double exponential fit in Equation (2.46). ??????????????????????????????.27 Fig. 2.6: Fante?s g(W) function plotted using the double exponential fit in Equation (2.46). ??????????????????????????????.28 Fig. 2.7: Normalized Correlation function using Equation (2.34) and Fante?s strong turbulence approximation of the Correlation function Equation (2.42). ?????29 Fig. 2.8: Relative scale sizes vs. propagation distance for an infinite plane wave. The point of intersection denotes the onset of strong fluctuations [8]. ???????.31 Fig. 2.9: Weak and Strong Variances compared with the Rytov value [18]. Circles shown represent the path of a composite curve that is predicted to be valid over all turbulence strengths. ????????????????????????...32 Fig. 2.10: Aperture averaging approximation using Fante?s weak correlation functions plotted for several variances or turbulence levels. ?????????????..37 Fig. 2.11: Aperture Averaging Factor using Fante?s strong correlation functions plotted for various C n 2 . ????????????????????????38 Fig. 2.12: Comparison of the Aperture averaging Factor using Fante?s weak and strong turbulence theory for a choice of variance = 4.32 which is considered as strong turbulence level. ??????????????????????????..39 Fig. 2.13: Comparison of the Aperture Averaging Factor using Fante?s weak correlation functions with a variance of 4.32e-3 implying weak turbulence and a viii second curve using Fante?s strong correlation functions with a variance of 4.32 implying strong turbulence. ?????????????????????...40 Fig. 3.1: University of Maryland Map showing the 863 meters link range used with the Transmitter shown on the roof of the AV. Williams Building and the receiver at the Chesapeake Building. ???????????????????????44 Fig. 3.2: Transmitter of the 863m link from AVW to Chesapeake Building. ???45 Fig. 3.3: Receiver of the 863m link from AVW to Chesapeake Building. ????.46 Fig. 3.4: Aperture Averaging Setup using a Schmidt-Cassegrain telescope with a 406.4 mm outer aperture size and a 127 mm inner obstruction diameter. The incoming beam intensity distribution is directed to a CCD camera to measure the variance of the received irradiance. ???????????????????46 Fig. 3.5: Captured frame using a shutter speed of 1/500 second showing turbulence scintillation effects. ?????????????????????????.47 Fig. 3.6: Frame analysis of received scintillation at the CCD camera. The dashed circles are different diameter regions selected to contain families of rays that enter magnified corresponding circular regions in the telescope aperture. ??????.49 Fig. 3.7: Front panel of the LabVIEW program designed to calculate and record irradiance statistics for the aperture averaging experiment. ?????????...50 Fig. 3.8: Experimental Data plotted versus the Andrews and Churnside Weak Turbulence Approximations (Equations (2.55) and (2.54) respectively). ????..53 Fig. 3.9: Aperture Averaging factor for Experimental Data with a variance of 0.44567 plotted along with Fante?s weak turbulence theoretical curve (using Equations (2.32), (2.33), (2.34)). ???????????????????????????.55 Fig. 3.10: Aperture Averaging factor for Experimental Data with a variance of 0.44567 plotted along with Andrew?s and Churnside weak turbulence approximations (Equation (2.55) and (2.54) respectively). ????????????????..56 Fig. 3.11 Relative Intensity Fluctuations (Aperture Averaging Factor) as a function of aperture diameter for varying . The points are experimental data for of 1.38 (Intermediate Turbulence) and the dashed lines and solid line are theoretical curves using weak and strong turbulence theory, respectively. ???????????.58 2 I ? 2 I ? Fig. 3.12: Experimental results plotted against derived Fante?s weak and strong turbulence theoretical curves for of 1.38. ???????????????.59 2 I ? ix Fig. 3.13: Aperture averaging factor for intermediate experimental data of point- detector variance 1.38 plotted along with the Andrews and Churnside weak turbulence approximations of the aperture averaging factor. ?????????.60 Fig. 3.14: Using Churnside Asymptotic theory for a spherical wave for all turbulence strengths (Equation (2.55)) to calculate the Aperture Averaging factor. The plane wave equation was used for the coherence length ? o with C n 2 = 3.1894e-14. ???61 Fig. 3.15: Frame Analysis using Rectangular apertures of the same area as the circular apertures chosen in the analysis. ????????????????????.63 Fig. 3.16: Aperture Averaging factor ?F? versus Aperture Area for Circular as well as rectangular apertures plotted for a variance of 0.44576. ???????????64 Fig. 4.1: Probability of detection and false alarm [8]. ????????????70 Fig. 4.2: Average BER versus S/N ratio in the absence of turbulence. Threshold is set to 1/2, and different curves are plotted for several intensity variances ( ). ???73 2 I ? Fig. 4.3: Effective S/N Ratio as a function of Aperture Diameter for an experimental variance of 1.38 plotted along with several theoretical variances calculated using Fante?s Correlation functions. ?????????????????????.77 Fig. 4.4: Effective signal-to-noise ratio for various aperture-averaged variances. ?.78 Fig. 5.1: Spherical/Bubble Model - Beam Wander Simulation. ????????.88 Fig. 5.2: Three dimensional Snell?s Law. ????????????????...92 Fig. 5.3: Flowchart of the Beam Wander Simulation Model. ?????????.94 Fig. 5.4: Mean Square Beam Wander and Cubic Fit for a 1km range using the Spherical/Bubble Model. ???????????????????????95 Fig. 5.5: Mean Square Beam Wander and Cubic Fit for a 5km range using the Spherical/Bubble Model. ???????????????????????96 Fig. 5.6: Spherical Bubble Model ? Phase Wander Simulation. ????????97 Fig. 5.7: Flowchart of the Phase Wander Simulation Model. ????????...100 Figure 5.8: Phase Wander Simulation Structure Function for Various Mean Free Paths. ??????????????????????????????.101 Fig. 5.9: Aperture Averaging Simulation using the Spherical Bubble Model. ??.103 x Fig. 5.10: Position Figure showing the origins of the rays in the x-y plane. ???104 Fig. 5.11a: Rays in the x-z plane. ???????????????????.106 Fig. 5.11.b: Rays in the y-z plane. ???????????????????106 Fig. 5.12: Random Three-dimensional Gaussian Ray propagation vector . ??..107 i r r Fig. 5.13: Ray Number verses Receiver Aperture Radius for 1 simulation run of Gaussian Distributed rays. Various curves are plotted for bubbles of varying standard deviation ?std?, n ? along with the no bubbles case. ????????????.111 Fig. 5.14: Ray Number verses Receiver Aperture Radius for 1 simulation run of Uniformly Distributed rays. Various curves are plotted for bubbles of varying standard deviation ?std?, n ? along with the no bubbles case. ???????......112 Fig. 5.15: Aperture averaging factor (F) versus the aperture radius for Gaussian distributed rays. Different curves are plotted for n ? =0.00001, 0.00005, and 0.0001. ?????????????????????????????..114 Fig. 5.16: Aperture Averaging Factor versus Aperture Radius for uniformly distributed rays. Different curves are plotted for several n ? values all within the strong turbulence regime. ??????????????????????..115 Fig. 5.17: Aperture Averaging Factor versus Aperture Radius for uniformly distributed rays. Different curves are plotted for several n ? values corresponding to different turbulence strengths. ????????????????????..116 xi List of Publications 1. Heba Yuksel, Stuart Milner, and Christopher C. Davis, "Aperture Averaging for Optimizing Receiver Design and System Performance on Free Space Optical Communication Links", Journal of Optical Networking, 4(8), 462-475, 2005. 2. Heba Yuksel, Walid Atia, and Christopher Davis, "A Geometrical Optics Approach for Modeling Atmospheric Turbulence," Paper presented at Conference 5891: Atmospheric Optical Modeling, Measurement, and Simulation, Optics & Photonics 2005, SPIE 50th Annual Meeting, San Diego, California, 31 July ? 4 August, 2005: To appear in the Proceedings of SPIE, Vol. 5891, 2005. 3. Heba Yuksel and Christopher Davis, "Aperture Averaging for Studies of Atmospheric Turbulence and Optimization of Free Space Optical Communication Links," Paper presented at Conference 5892: Free-Space Laser Communications V, Optics & Photonics 2005, SPIE 50th Annual Meeting, San Diego, California, 31 July ? 4 August, 2005: To appear in the Proceedings of SPIE, Vol. 5892, 2005. Conference Papers: 1. Heba Yuksel, Christopher C. Davis, Linda Wasiczko, "Aperture Averaging Experiment for Optimizing Receiver Design and Analyzing Turbulence on Free Space Optical Communication Links", Paper presented at CLEO 2005, Baltimore, Maryland, May 22-27, 2005. OSA Technical Digest (OSA, Washington DC 2005). 2. Heba Yuksel, Christopher C. Davis, "Optimizing Receiver Design Through Turbulence Analysis of Free Space Optical Communication Links", Paper presented at the Optics in the SouthEast Conference'2004, Conference Technical Digest p. 91, Charlotte, North Carolina. xii Chapter 1: Introduction 1.1 Thesis Contributions This dissertation shows a flexible empirical approach for improving link performance through image analysis of intensity scintillation patterns coupled with video-frame aperture averaging on a free space optical (FSO) communication link. An imaging system for measuring the effects of atmospheric turbulence and obscuration on FSO links will be presented. Weak and intermediate turbulence results will be shown for an 863 meter link at the University of Maryland. In addition, efficient computational techniques have been developed for various correlation functions that are important in assessing the effects of turbulence. This thesis will present the most accurate empirical data to date for the intermediate turbulence regime. Such results can help develop upon existing empirical data and lead to the development of new theories. In addition, the weak turbulence results show the best fit to date to the analytical expressions described by theory. A new geometrical model to assess the effects of turbulence on laser beam propagation will also be presented. The atmosphere along the laser beam propagation path is modeled as a spatial distribution of spherical bubbles with refractive index discontinuities that are statistically distributed in size and diameter according to various models. The Monte Carlo techniques used allow us to assess beam wander and phase shifts effects along the path, and aperture averaging effects at the receiver. The Monte Carlo simulations are compared and are well fitted with the predictions of 1 wave theory. This thesis is the first to date to geometrically simulate aperture averaging. 1.2 Turbulence Overview Direct line-of-sight optical communication links, which are commonly called ?optical wireless? systems or free space optical (FSO) communication links are becoming increasingly popular. Such links can provide virtually unlimited bandwidth at a relatively low cost and high performance communication over short distances up to a few kilometers. In addition, FSO links do not require any spectrum allocation by the Federal Communications Commission (FCC) and are highly secure because of their directionality with a low probability of interception and/or detection. Such links are also rapidly deployable, scalable and flexible. These properties make them attractive in a variety of applications including communication between buildings in cities and industrial parks, and between aircraft or ships, especially for military purposes due to their high security of transmission. Unfortunately, the atmosphere is not an ideal communication channel. Atmospheric turbulence can cause fluctuations in the received signal level, which increase the bit errors in a digital communication link. In order to quantify the performance limitations, a better understanding of the effect of the intensity fluctuations on the received signal at all turbulence levels is needed. The local density of the atmosphere is constantly fluctuating because of temperature and pressure fluctuations. This is atmospheric turbulence. The foundations of the study of atmospheric turbulence were laid in the late 1960s and 1970s. Theory reliably describes the behavior in the weak turbulence regime, but 2 theoretical descriptions in the intermediate and strong turbulence regimes are less well developed. Within the troposphere layer of the atmosphere, the air temperature decreases rapidly with increasing altitude causing turbulence strength to decrease with altitude. Turbulence can become very strong near the ground because of the heat transfer between ground and air. However, even at a modest height of 12 meters, which is used in the experiment to be presented, strong turbulence is rarely experienced. Nevertheless, strong turbulence can still occur if the beam propagates near rooftops and/or near other structures that cause increased temperature fluctuations. In addition, strong turbulence effects can occur for a beam leaving and entering an aircraft because of boundary-layer turbulence. An FSO system is expected to perform in both weak and strong turbulence conditions. Because of the shortcomings of theory, an empirical approach to determining ways of improving link performance through experiment and data analysis will fill gaps in the theory, and potentially result in the development of new theories. When a laser beam propagates through the atmosphere the randomly varying spatial distribution of refractive index that it encounters causes a number of effects. These include: (1) A fluctuating intensity as observed with an optical detector at the end of the path. This is referred to as scintillation. (2) A varying degree of fluctuation with the size of the detector, or with the size of the receiving optics that direct the collected light to the detector. This is called aperture averaging. 3 (3) If a circularly symmetric Gaussian beam is observed at different distances from a transmitter it suffers progressive deterioration with increasing distance and turbulence strength. The progressive changes that are observed are: (i) deviations of the beam shape from circular that are time dependent; (ii) wander of the centroid of the beam; (iii) increase in the width of the beam over and above that expected from diffraction; and (iv) breakup of the beam into distinct patches of illumination whose shapes and locations fluctuate with time. (4) The ?coherence length? of the laser beam falls. (5) The angle of arrival of the phase fronts at a receiver fluctuates. Winds, which move the atmosphere in a more correlated way, can cause the centroid of the beam to shift, but they do not intrinsically randomize the laser beam as does turbulence. In principle, the effects of correlated atmospheric effects that ?steer? a laser beam can be compensated with a beam-steering scheme at the transmitter. 1.3 Aperture Averaging Intensity fluctuations at a receiver lead to a received power variance that depends on the size of the receiver aperture. Increasing the size of the receiver aperture reduces the power variance. This effect of the receiver size on power variance is called aperture averaging and will be presented and studied in detail in this dissertation. If there were no aperture size limitation at the receiver, then there would be no turbulence-induced scintillation. In practice, there is always a tradeoff between aperture size, transceiver weight, and potential transceiver agility for pointing, 4 acquisition and tracking (PAT) of links. An optimized receiver aperture size can be selected by using multi-frame image analysis of received intensity scintillation patterns. Aperture averaging theory has been extensively developed for plane and spherical waves in weak turbulence conditions [1-7]. Minimal theory is available for the strong turbulence regime [6, 8]. There has been some previous experimental work [2, 3], but early experiments did not account for scintillation saturation, and resulted in data that is different from that predicted by theory. Later data were limited by the short path lengths under investigation [3]. The Maryland Optics Group (MOG) at the University of Maryland has been developing novel technologies that will enable autonomous reconfiguration in hybrid, wireless networks with directional optical and RF links [9]. These include: 1) autonomous connection/configuration of narrow beam links; and 2) management of changes in link states in the context of an overall network. The former is referred to as the pointing, acquisition and tracking (PAT) of links, and the latter as topology and link control. In this thesis, a flexible, empirical approach is described for optimizing the design and configuration of FSO communication links by using multi-frame image analysis of received intensity scintillation patterns. This is a versatile way for performing aperture averaging analysis, which guides receiver aperture size selection when there are trade offs between transceiver size, weight, and power. The Maryland Optics Group has previously described the use of hybrid FSO/RF systems to provide better directional wireless connectivity than either FSO or RF alone [9]. In clear conditions, FSO links can provide significant power margins for excellent performance, even over long distances, but these links must be able to 5 operate in the face of potentially large fluctuating received powers. These ?fades? and ?surges? result from the imperfect nature of the optical communication channel through the atmosphere caused by atmospheric turbulence. Atmospheric turbulence is the random fluctuations of coupled temperature/density/refractive index, which vary from point to point along the link path. A key to the maintenance of network quality of service, amidst such degradation, is to optimize the level of received signals [10]. In this thesis, aperture averaging of received scintillation will be studied in order to design receivers that can mitigate the negative effects of turbulence and hence optimize performance. Intensity scintillation and associated beam wander are the two most problematic characteristics of an FSO link. A receiver must be large enough to collect sufficient power and reduce scintillation effects at a given range, but must also be small enough to be of practical size for cost effectiveness. In general, an FSO receiver does not need to have a very large diameter to reduce scintillation effectively, because of a nonlinear reduction in scintillation with increasing the aperture size. Optimized apertures increase the effective signal-to-noise (S/N) ratio and reduce the bit-error-rate (BER). 1.4 Geometrical Simulations In line-of-sight optical communication systems, absorption of the beam by the atmosphere can be important, especially in adverse weather conditions of fog, snow, heavy rain, or in conditions of battlefield obscuration. The combined effects of direct absorption and scattering of laser light can be described by a single path-dependent absorption coefficient ()z? . The power reaching a receiver (RX) from a transmitter (TX) is easily calculated for links without significant turbulence effects. The received 6 power for a RX area A, range L, and beam divergence angle ? varies as 22 ? ? LPAe L? ? , where P is the TX power, and ? in this case is assumed constant along the path L. When turbulence effects are included, the effects of the atmosphere are in a sense more subtle. The laser beam diverges as it propagates, which leads to a reduction in the received power with range. In addition, the phase fronts of the laser beam are distorted, and in coherent optical communication applications, in which the received beam is mixed with a local oscillator laser beam -- much as in a FM radio -- badly distorted phase fronts lead to poor mixing and serious loss of detectable signal. In line-of-sight optical links in which the receiver acts as a ?photon-bucket? for received light -- only the amount of received optical power is important, not the quality of the incoming wavefronts. However, indirectly, the effects of the atmosphere on phase fronts does lead to ?beam wander? and problems associated with pointing the transmitter beam at the receiver, so this issue can not be completely ignored. In principle, such effects can be minimized by ?steering? the beam so that it always hits the receiver. Even when this is done, however, several problems remain, even for a photon-bucket or intensity-based receiver. This is particularly true for a terrestrial link that runs essentially parallel to the ground. Links from ground to space, and vice-versa, are affected differently because the laser beam travels for only a small part of its path through the denser layers of air near the ground where atmospheric disturbances are concentrated. In this thesis, a theoretical model is presented with new features for calculating the effect of atmospheric turbulence and obscuration on line-of-sight laser communication links. This model will build on existing weak and strong turbulence 7 theories to assess the effect of the turbulence parameter on link performance where z is distance along the propagation path. Efficient computational techniques are presented in this dissertation for various important correlation functions that are important in assessing the effects of turbulence. An additional, new geometrical optics approach is also developed in which the atmosphere is modeled as a spatial distribution of spherical bubbles in which there are small refractive index distributions distributed according to various statistical models. The model has proved capable of assessing beam wander, phase shifts, and aperture averaging at the receiver for a laser beam propagating down range in turbulence. A receiver of varying size can be used to collect rays that are initially Gaussian distributed or uniformly distributed and propagate through the simulated spherical bubbles. This allows aperture averaging to be calculated for varying receiver sizes. The random nature of the atmosphere along the link is constructed from random number generated distributions, and the effective can be determined for a given randomly selected atmospheric path by correlating beam wander behavior with the path length. These statistical analyses are obtained by Monte Carlo simulations over distributions of randomly selected index distributions along atmospheric paths. )( 2 zC n 2 n C 1.5 Thesis Organization This thesis is divided into 6 chapters. Chapters 1 and 6 are the introduction and conclusion respectively. Chapter 2 presents an overview of the key turbulence theory that will be used throughout this thesis. It also includes efficient computational techniques for various important correlation functions that are important in assessing 8 the effects of turbulence. Chapter 3 presents the experimental setup, methodology, and results of an imaging system that is used to measure the effects of atmospheric turbulence and obscuration on FSO links. Results are presented for weak and intermediate turbulence and compared with the atmospheric turbulence theory developed and presented in Chapter 2. Chapter 4 describes the aperture averaging effect on the performance of FSO links through bit-error-rate (BER) and Signal-to- noise (S/N) ratio analysis. Finally, Chapter 5 describes a new geometrical optics model than can assess the effect of atmospheric turbulence on the propagation of rays that are initially Gaussian as well as uniformly distributed. 9 Chapter 2: Turbulence Theory Overview 2.1 Introduction A complete theory of wave propagation through random media is not yet available ? it remains an active area of research in many diverse fields such as atmospheric optics, ocean acoustics, radio physics, and astronomy [4]. However, the general theory is fairly well understood in certain asymptotic regimes. We will concentrate in this chapter on a plane wave approximation for the propagating beam. 2.2 Atmospheric Turbulence Optical Turbulence can be defined as the fluctuations in the index of refraction resulting from small temperature fluctuations. The atmosphere has two distinct states of motion as a viscous fluid ? laminar and turbulent. In the earliest study of turbulent flow, Reynolds defined a non-dimensional quantity Re=Vl/v, where V is the characteristic velocity (in m/s), l is the dimension of the flow (in m), and v is the kinematic viscosity (in m/s 2 ). When the flow of a viscous fluid exceeds a critical Reynolds number, the flow changes from laminar to a more chaotic state called turbulence. Turbulent air motion is represented by a set of eddies of various scale sizes extending from a large scale size L o called the outer scale of turbulence to a small scale size l o called the inner scale of turbulence forming the inertial range. Scale sizes smaller than the inner scale belong to the dissipation range. The source of 10 energy at large scales is either wind shear or convection. Under the cascade theory, the wind velocity increases until it reaches a point at which the critical Reynolds number is exceeded. This causes local unstable air masses called eddies to form. Under the influence of inertial forces, larger eddies break up into smaller eddies forming the inertial range between scale size of l o to L o . Figure 2.1 shows the Kolmogorov cascade model of turbulence as a function of spatial scale [4]. Turbulence Spectrum ?(?) Energy Input Energy Dissipation ? 0 Wavenumber ? ? s Fig. 2.1: A pictorial description of the process of turbulent decay. As turbulent eddies subdivide, they become smaller and more uniform until all of their energy dissipates as heat [4]. The outer scale L o is usually assumed to grow linearly with the order of the height above the ground of the observation point up to approximately 100 meters. Eddies of scale sizes smaller than L o are assumed statistically homogeneous and isotropic. Statistical homogeneity implies that the mean value of the field is constant and that correlations between random fluctuations in the field from point-to-point are independent of the chosen observation points, depending only on their vector separation. Statistical isotropy implies that point-to-point correlations depend only on 11 the magnitude of the vector separation between observation points [4]. The inner scale l o is typically on the order of 1 to 10mm near the ground but is in the order of centimeters or more in the troposphere and stratosphere. The most important parameter to consider in optical wave propagation is the index of refraction fluctuations caused by atmospheric turbulence [4], )()( 1 rnnrn o rr += , (2.1) where r r is a point in space, 1)( ?= rnn o r is the mean value of the index of refraction of air at atmospheric pressure and )( 1 rn r represents the random deviation of )(rn r from its mean value. Therefore 0)( 1 =rn r . Time variations are suppressed in the calculations assuming the wave maintains a single frequency as it propagates. The fluctuations in the index of refraction are related to corresponding temperature and pressure fluctuations as follows [4], )( )( 10791)( 6 rT rP rn r r r ? ?+= , (2.2) where P is the pressure in millibars and T is the temperature in degree kelvin. Humidity fluctuations only contribute in the far-IR region, and pressure fluctuations are usually negligible. Therefore, the index of refraction fluctuations within the visible and near-IR region of the spectrum are due primarily to random temperature fluctuations [4]. The covariance function of )(rn r can be expressed by, 2 11111121 )()(),(),( onn nrrnrnrrrBrrB ++=+= rrrrrrrr , (2.3) 12 where and are 2 points in space and 1 r r 2 r r 12 rrr rrr ?= . Assuming a homogeneous, isotropic turbulent media, the covariance function reduces to a function of only the scalar distance 12 rrr rr ?= . Random fields that permit a decomposition into a varying mean and a statistically homogeneous fluctuation are called locally homogeneous [4]. Locally homogeneous fields are usually not characterized by the covariance function, but by the structure function, ()( )[][ ])()0(2rrr)( 2 11 rBBnnrD nnn ?=?+= . (2.4) By plugging Equation (2.3) into Equation (2.4), we arrive at the Kolmogorov- Obhukov two-thirds power law describing the structure function of refractive index fluctuations for separations between l 0 and L 0 [4], 322 )( rCrD nn = , oo Lrl <<<< (2.5) where is the index of refraction structure parameter, also called the structure constant. is a function of height. Over short time intervals at a fixed propagation distance and constant height above the ground which is the case of interest in this research, it is reasonable to assume that is essentially constant. typically ranges from 2 n C 2 n C 2 n C 2 n C 3217 10 ?? m or less for conditions of ?weak turbulence? and up to 3213 10 ?? m or more when the turbulence is ?strong?. 2.3 Power Spectrum Models for Refractive index fluctuations We have described thus far the covariance as well as the structure function of the index of refraction fluctuations. The three-dimensional spatial power spectrum of 13 the random field )(? r n ? form a Fourier transform pair with the covariance function [11], , (2.6) ??? ? ?? ?= ?? ? 3. )()( derB n ri n rr rr ??? ? ?? ? ? ? ? ? ? ? =? rdrBe n ri n 3. 3 )( 2 1 )( rr rr ? ? ? , (2.7) where ? r is the wave number vector. Assuming homogeneity and isotropy, these Fourier transform relations reduce to [4], ? ? =? 0 2 )sin()( 2 1 )( rdrrrB nn ? ?? ? r , (2.8) ? ? ?= 0 )sin()( 4 )( ???? ? dr r rB nn r , (2.9) where l/2??? == r is the magnitude of the wave number vector and l is the turbulent eddy size. If we call the one-dimensional spectrum )(? n V , then its relation with the three-dimensional spectrum )(? r n ? is given by [4], ? ? ?? ? d dV n n )( 2 1 )( ?=? . (2.10) Therefore a one-dimensional spectrum exhibiting a 35? ? behavior corresponds to a three-dimensional spectrum with a 311? ? behavior. In addition, the relation between the structure function and the power spectrum is given by [4], 14 ? ? ? ??? d r r rD nn ? ? ? ? ? ? ? ? ??= 0 2 )sin( 1)(8)( , (2.11) drrD dr d r dr d r r nn ? ? ? ? ? ? ? ? =? 0 2 22 )( )sin( 4 1 )( ? ? ?? ? . (2.12) The correlation and covariance functions represent a spatial domain description whereas the power spectrum is a wave number representation. 2.3.1 Kolmogorov Spectrum Kolmogorov defined the power spectral density for refractive index fluctuations over the inertial range by [4], , 3/112 033.0)( ? =? ?? nn C 0 /1/1 lL o ??? . (2.13) The Kolmogorov spectrum is used when the inner scale is zero and the outer scale is infinite, or as long as the wave number is within the inertial subrange (where eddysize/2?? = ). Turbulence Spectrum ?(?) Inertial subrange Fig. 2.2: Spectrum of the refractive index fluctuation. The energy input range, inertial subrange, and energy dissipation ranges are indicated [11]. ? 0 =2?/L 0 Wavenumber ? (~? -11/3 ) Energy Energy dissipation Input range range ? 0 =2?/l 0 15 2.3.2 Tatarski Spectrum When the inner or outer scale effect cannot be ignored, the Kolmogorov power law spectrum in Equation (2.13) needs to be modified. Tatarski suggested the extension of Equation (2.13) into the dissipation range o l/1>? through the introduction of a Gaussian function that essentially truncates the spectrum at high wave number [12], , )/exp(033.0)( 223/112 mnn C ???? ?=? ? o L/1?? , (2.14) where om l/92.5=? . This method was developed for the purpose of mathematical convenience. However, Equation (2.14) suggests a singularity at 0=? for the limiting case . This implies that the structure function (Equation 2.11) exists but the covariance function (Equation 2.9) does not. 0/1 = o L )(rD n )(rB n 2.3.3 Von Karman Spectrum Von Karman modified the Tatarski spectrum such that it is finite for o L/1 R I ? 2 ln R I ? A second important parameter is the variance of the normalized intensity fluctuations, which is [11], 2 2 2 2 I II I ? =? . (2.1) For weak turbulence the Rytov variance can be re-written as [18], 2 2 ln ln ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? = I I R I ? 2 1 ? ? ? ? ? ? ? ? ?= I I (2.2) 2 2 2 I II ? = . Consequently, in weak turbulence [18], . (2.3) 22 ln II R ?? = In other words, if the intensity variance is observed under conditions of weak turbulence, it will be identical to the variance of the log intensity. 2.4.2 Numerical Calculation of Fante?s Correlation Functions in Weak Turbulence Although theoretically there are differences in the way the atmosphere perturbs plane waves, spherical waves, and focused laser beams, there is considerable similarity between many of the effects on plane waves and collimated laser beams. 19 Consequently, only plane waves will be dealt with specifically in this discussion and in the case of laser beam wander, collimated laser beams. The field of such a wave propagating in the z direction can be represented at the transmitter as [18] )( 0 EE kztj e ? = ? . (2.4) The magnitude of the wave vector k is k = 2? /? , where ? is the wavelength of the wave. For a plane wave E 0 is constant over the transmitter aperture, which we will assume is located in the plane z = 0. For a linearly polarized wave E always points in the same direction. For a laser beam the field distribution at the transmitter is [18] 22 0 E)(E wr er ? = , (2.5) where w is the spot-size at the transmitter. Because the atmosphere is not intrinsically chiral, left and right circularly polarized waves should be identically affected by turbulence so we do not expect any perturbation of the polarization state of a light wave that has propagated through turbulence. At a receiver aperture located at z = L, the electric field fluctuates in time and space because of turbulence. Because turbulence is a random phenomenon, the actual behavior of the electric field E(z, r, t) components of a wave propagating through turbulence cannot be determined. All that can be calculated are various time and ensemble averages over field variables. Of central importance is the correlation function of log amplitudes for two observation points spaced a distance ? apart in a plane that is perpendicular to the direction of wave propagation, ( ) ( ) ( ) ,,,, 21 ????? ? zzzB = (2.26) 20 where 21 ??? ?= . The distance z is measured along the line of sight from the source. For weak turbulence, provided the log amplitude? is normally distributed, this correlation coefficient can be related to the normalized correlation for intensity ()?,zb I defined by the relation [4, 18], () ()( ) ( ) ( ) ()() , ,, ,,,, , 21 2121 ?? ???? ? zIzI zIzIzIzI zb I ? = (2.27) by () ()[?? ? ,1ln 4 1 , zbzB I += ]. (2.28) Calculations of these correlation coefficients differ somewhat depending on whether a plane wave or Gaussian beam is propagating from source to observation plane, and on whether the Gaussian beam is being focused or not. I will deal here only with the plane wave case: a large diameter Gaussian beam that is diverging slowly is equivalent to a plane wave provided w<>+= I R I L ? ? ? , (2.36) where is the Rytov variance. 2 R ? For the spherical wave case, the Andrews asymptotic analysis used the Kolmogorov spectrum to relate the irradiance variance in the saturation region to the plane wave Rytov variance [3], 1, )( 73.2 1)( 2 5/22 2 >>+= I R I L ? ? ? . (2.37) 2.5.2 Churnside Asymptotic Analysis After the phenomenon of saturation of scintillation was understood, Churnside built upon Fried?s work and published the first significant application of asymptotic theory to the study of aperture averaging. Churnside?s result for the irradiance variance in strong turbulence for the plane wave case is, 1, )( 14.1 1)( 2 5/22 2 >>+= I R I L ? ? ? . (2.38) For the spherical wave case, the Churnside approximate of the scintillation index in strong turbulence is, 3/1 2 2 86.31 ? ? ? ? ? ? ? ? += L k o I ? ? , (2.39) where the coherence length for the spherical wave o ? is, 5/322 )546.0( ? = no LCk? . (2.40) The spherical wave scintillation index in strong turbulence can be expressed in terms of the plane wave Rytov variance as, 25 1, )( 35.5 1)( 2 5/22 2 >>+= I R I L ? ? ? . (2.41) 2.5.3 Numerical Calculation of Fante?s Correlation Functions in Strong Turbulence In strong turbulence, for which , the correlation function for normalized intensity is [19], 1 2 ln >> R I ? () () () () , 1 2exp, 2/1 11/6 ln 5/6 ln 2/15/2 2 ln 3/5 0 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? ? ? ? ? ? ? ? + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?== z g z fLzb R R R I I I I ? ?? ?? ? ? ? ? ? (2.42) where is the plane wave lateral coherence length for uniform turbulence along a path of length L. The correlation function of log amplitudes can be calculated from Equation (2.34). The functions ( 5/3 22 0 46.1 ? = n LCk? ) ( )Sf and ( )Wg used in Equation (2.42) in evaluating the correlation function for normalized intensity in strong turbulence have been given by Fante [23, 24], although beware of a critical error in the power of t in the g(W) function in Equation (6) of reference [24], where it is written as t 8/3 but should be t 8/11 as defined correctly by Fante in reference [23]. These functions are () () ( ?? ? ???? = 1 00 15/3 0 66.226.45/23/1 ,54.343.1 SytJedttdyySf yt ) ) (2.43) and () ( ) ( ?? ?? ??? ?= 00 11/311/8 0 3/8 .43.2cos127.0 WstJdsetdttWg s (2.44) 26 The function is easier to evaluate numerically than()Sf ( )Wg and the range of integration does not need to be split into as many parts as it is for in order to obtain a convergent answer. For example, the inner integral in Equation (2.43) is reliably evaluated by breaking the integration range into parts from 0 to 1 and 1 to 10. On the other hand the double integral in Equation (2.44) was broken up into 12 subranges to obtain a reliable result. For use in calculating correlation functions, such as Equation (2.42) it is convenient to have a simpler version of ()Wg ( )Sf and . The function is ?very well? fit by a tenth order polynomial in the range with ()Wg ()Sf 50 ?? S () 1098 7654 32 0000220974.00006918.00093015.0 0701426.03241827.0937165.06331602.1 4652269.11696313.07527889.04203171.0 SSS SSSS SSSSf ?+? +?+? +??= . (2.45) -0.05 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 0.45 0123456 S f(S ) Fig. 2.5: Fante?s f(S) function plotted using the tenth order polynomial in Equation (2.45). The function is very well fit by a double exponential function in the range with ()Wg 100 ??W 27 () .3786.01852.0 2818.10744.0 WW eeWg ?? += (2.46) 0 0.1 0.2 0.3 0.4 0.5 0.6 02468101 W g(W ) 2 ? Fig. 2.6: Fante?s g(W) function plotted using the double exponential fit in Equation (2.46). A check on the validity of the numerical results is provided from the value , which is numerically determined to be 0.5633. An analytic result for can be determined since ()0g ()0g () ( ) ? ?? ?? ?= 00 3/8 ,cos127.00 dsetdttg s (2.47) which gives () ( ) ? ? ? ?= 0 3/8 .cos127.00 dtttg (2.48) Integration by parts and using the result from Gradshteyn and Ryzhik [25] that ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?= 0 3/2 , 6 cos 3 1 cos ? tdtt (2.49) 28 gives , which is in excellent agreement with the numerical result. () 563767.00 =g Figure 2.7 plots the normalized correlation function (Equation (2.34)) using strong turbulence theory (Equation (2.42)), our link parameters and the f(S) and g(W) fitted functions in Equations (2.45) and (2.46) respectively. 0 0.5 1 1.5 2 2.5 0 0.05 0.1 0.15 0.2 0.25 ? (m) N o r m a l iz e d C o r r e l a t ion F unc t i on, b I ( ? ) ( S t r ong The or y ) Cn2=1e-12 Cn2=1e-13 Cn2=1e-14 Fig. 2.7: Normalized Correlation function using Equation (2.34) and Fante?s strong turbulence approximation of the Correlation function Equation (2.42). The aperture averaging measurements in this thesis will be compared with the predictions of weak and strong turbulence theory. 2.6 Weak to Strong Turbulence Theory As a coherent wave propagates in the atmosphere, the wave is scattered by the smallest of the turbulent cells (on the order of millimeters) through diffraction. The largest turbulent cells within the inertial range act as refractive ?lenses? with focal lengths typically on the order of hundreds of meters or more. Small-scale 29 contributions to scintillation are associated with turbulent cells smaller than either the first Fresnel zone kL / or the coherence radius o ? , whichever is smallest. The Fresnel zone defines the most effective turbulent cell size in producing scintillation at distance L from the source. Turbulent cell sizes smaller than the Fresnel zone contribute less to scintillation because of the weaker refractivity fluctuations associated with them, and cell sizes larger than the Fresnel zone do not diffract light through a large enough angle to reach the receiver at L. Large-scale fluctuations of the irradiance are generated by turbulent cells larger than that of either the Fresnel zone or the scattering disk o kL ?/ , whichever is largest. The scattering disk is defined by the refractive cell size l at which the focusing angle Ll F /?? is equal to the average diffraction angle oD k?? /1? [8]. Figure 2.8 shows the relative scale sizes versus the propagation distance L for an infinite plane wave, ( ) 5/3 22 46.1 ? = LkC no ? with wavelength m?? 06.1= , a fixed , and inner-scale and outer-scale effects ignored. The onset of strong fluctuations occurs just beyond 200 m where the curves intersect. For constant and weak fluctuations, the scale size of the spatial coherence radius 3/2132 105 ?? = mxC n 2 n C o ? is larger than the Fresnel zone size, but the Fresnel zone represents the correlation width of the optical wave and is the most effective cell size in producing irradiance fluctuations in this regime. At the onset of moderate-to-strong fluctuations, the spatial coherence radius approaches the scale size of the Fresnel zone, and hence all three cell sizes (spatial coherence radius, Fresnel zone size, and scattering disk) are roughly equal. This happens in the vicinity of the focusing regime where irradiance fluctuations are 30 maximum. For stronger irradiance fluctuations, the spatial coherence radius and the scattering disk are the dominant cell sizes forming the upper bound of small diffractive cells and lower bound of the large refractive cells, respectively. Cell sizes between that of the coherence radius and the scattering disk have little effect on scintillation in this regime [8]. 0.0 0.5 1.0 1.5 2.0 2.5 0 100 200 300 400 500 Distance (m) S c al e S i z e ( c m) Coherence Radius Scattering Disk Fresnel Zone ? =1.06 ? m C n 2 = 5x10 -13 m -2/3 Fig. 2.8: Relative scale sizes vs. propagation distance for an infinite plane wave. The point of intersection denotes the onset of strong fluctuations [8]. Experimental data reveals that the scintillation index increases initially within the weak turbulence regime with increasing values of the Rytov variance. It then reaches its maximum value in the focusing regime and gradually decreases toward a unity as the Rytov variance increases without bound. 31 The value of , which is the autocorrelation function of the intensity observed at a particular point, is the same as the intensity variance that will be observed with a point detector. Figure 2.9 compares the variances calculated using Fante?s weak and strong turbulence theories with the Rytov value [18]. It also shows the path of a composite curve that is predicted to be valid over all turbulence strengths. )0,(Lb I COMPARISON OF CORRELATION PARAMETER WITH RYTOV VALUE Square root of variances, 863 meters range, 632.8 nonometers 0 1 2 3 012345 Sigma(Rytov) S i g m a( W e ak ) and S i gm a ( S t r o n g ) Rytov Weak Theory Strong Theory Composite Weak Rytov Strong Fig. 2.9: Weak and Strong Variances compared with the Rytov value [18]. Circles shown represent the path of a composite curve that is predicted to be valid over all turbulence strengths. 32 Note that the ?Weak? curve follows the linear variation predicted by the Rytov variance well up to . The ?Strong? curve is not valid for small values of , but should represent expected variances well for . The circles in Figure 2.9 represent a path for a ?composite curve? that follows 3.0 2 ln = R I ? 2 ln R I ? 1 2 ln > R I ? weak ? for small and then merges with 2 ln R I ? strong ? for larger . Such a composite curve demonstrates the increase in intensity variance that occurs as turbulence increases followed by its saturation at high levels of turbulence. Ishimaru [11], Fig.(20-10) presents such a ?schematic curve?, and shows the saturation of the intensity variance at a value just above unity. 2 ln R I ? 2.6.1 Scintillation Index Model Andrews and Phillips developed a model for the plane wave case that is valid under all fluctuation conditions. This was done through replacing the Kolmogorov spectrum (Equation (2.13)) with the effective Kolmogorov spectrum [7, 8, 26], [ ] [)()(033.0 )()()()( 3/112 , ??? ] ???? yxn yxnen GGC GG += +?=? ? , (2.50) where the large-scale filter function that passes only spatial frequencies x ?? < is, ? ? ? ? ? ? ? ? ?= 2 2 exp)( x x G ? ? ? , (2.51) and the small-scale filter function passing only y ?? > is, 6/1122 3/11 )( )( y y G ?? ? ? + = . (2.52) 33 Using the modified Rytov theory, the scintillation index is defined in terms of the large-scale and small scale log irradiance fluctuations which are described in detail in Reference [7], 2 ln x ? 2 ln y ? . (2.53) 1)exp( 2 ln 2 ln 2 ?+= yxI ??? The scintillation index for a plane wave, excluding inner scale effects, is then defined to be [7, 8], 1 )69.01( 51.0 )11.11( 49.0 exp 6/55/12 2 6/75/12 2 2 ? ? ? ? ? ? ? + + + = R R R R I ? ? ? ? ? , (2.54) where is the Rytov variance for a plane wave. Equation (2.54) reduces to the Rytov variance under weak turbulence conditions, and to the Andrews asymptotic model in Equation (2.36) under strong turbulence conditions. 2 R ? For a spherical wave, the scintillation index for zero inner scale effects is defined by [7, 8], 1 )23.01( 20.0 )19.01( 20.0 exp 6/55/12 2 6/75/12 2 2 ? ? ? ? ? ? ? + + + = R R R R I ? ? ? ? ? . (2.55) 2.7 Aperture Averaging Intensity fluctuations at a receiver lead to a received power variance that depends on the size of the receiver aperture. Increasing the size of the receiver aperture reduces the power variance. This effect of the receiver size on power variance is called aperture averaging. The aperture averaging factor ?F? is defined as the ratio of the normalized intensity variance of the signal at a receiver with diameter D to that of a point receiver [27], 34 . )0( )( 2 2 = = D D F I I ? ? (2.56) Received signals are normalized by the square of the average signal. Tatarski has given an expression that allows the aperture averaging effect to be calculated from the correlation function of the normalized intensity [27]. His result is ( ) () () ? ? = 0 2 , 0 16 ??? ? ? dK b b D F I I (2.57) where b I (?) is the covariance function of the irradiance, b I (0) is the variance of the irradiance, which is equivalent to the Rytov variance for weak turbulence, D is the diameter of the receiver aperture, and () ( ) ( ) ( )[ ] 21 22 1arccos DDDK ???? ??= . (2.58) The aperture averaging factor F represents the intensity variation seen with the actual receiver relative to a point receiver. In practice, at range L, a receiver whose diameter satisfies LD ?<< will behave as a point receiver. The calculation of F requires an integration involving the correlation function ( )? I b and it is in this calculation that the simple expressions for this function in weak and strong turbulence (Equations (2.33), (2.34), (2.42), (2.45), (2.46)) become valuable. For a plane wave with small inner scale, l 0 << (L/k) 1/2 , Churnside approximated the aperture averaging factor in weak turbulence by [3, 6, 7], . 4 07.11 1 6/7 2 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? += L kD F (2.59) In weak turbulence, Andrews reported a better approximation for the aperture averaging factor as [4], 35 . 4 062.11 6/7 2 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? += L kD F (2.60) For small apertures (kD 2 /4L << 1), F=1 as expected. For larger apertures (kD 2 /4L >> 1), the variance decreases with increasing aperture size. Note that as turbulence becomes more severe, aperture averaging becomes more effective in reducing the intensity variance, but only up to a point. Significant aperture averaging kicks in at a very small receiver diameter, but there is a long tail. This effect occurs because strong turbulence scrambles the beam sufficiently that it becomes almost homogeneous and reduces the intensity variance. In strong turbulence, Churnside approximated the strong turbulence spherical wave aperture averaging factor as [3], 1 3/7 2 2 1 2 2 2 2 613.01 2 1 2 908.01 2 1 ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + + = L kD D F o I I oI I ? ? ? ?? ? , (2.61) where is the scintillation index in strong turbulence conditions defined by Churnside in Equations (2.38), (2.41) and Andrews in Equation (2.36), (2.37) and the transverse coherence length for the spherical wave is . 2 I ? 5/322 )545.0( ? = no LCk? 2.7.1 Aperture Averaging calculations using Fante?s correlation functions in weak turbulence Using our calculations of the normalized correlation functions in weak turbulence (Equations (2.33) and (2.34)), the aperture averaging factor can be calculated using Equation (2.57). Figure 2.10 plots the approximation of the aperture 36 averaging factor in weak turbulence for several variances. Only the variance of value 4.23e-3 is considered a weak turbulence variance. The rest of the variance values plotted actually correspond to intermediate to strong turbulence variances and such a theoretical equation used in plotting the curves applies only for weak turbulence regimes and will provide inaccurate measurements for non-weak variance curves. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 051015 Aperture Diameter (in cm) A p er t u r e A v er ag i n g F a ct o r , F 1e-13 (4.32) 5e-14 (2.16) 1e-14 (0.432) 1e-16 (4.32e-3) C n 2 (? I 2 ) Fig. 2.10: Aperture averaging approximation using Fante?s weak correlation functions plotted for several variances or turbulence levels. 37 2.7.2 Aperture Averaging calculations using Fante?s correlation function in strong turbulence Using our calculations of the normalized correlation functions in strong turbulence (Equations (2.34) and (2.42)), the aperture averaging factor can be calculated using Equation (2.57). Figure 2.11 plots the approximation of the aperture averaging factor using strong turbulence theory for various values of C n 2 or turbulence strengths. 0 0.2 0.4 0.6 0.8 1 0 0.05 0.1 0.15 0.2 Diameter (in m) A p er t u r e A ver ag i n g F act o r , F 1e-12 1e-13 1e-14 C n 2 Fig. 2.11: Aperture Averaging Factor using Fante?s strong correlation functions plotted for various C n 2 . Figure 2.12 plots the approximation of the aperture averaging factor using weak and strong turbulence theory for the same C n 2 of 1e-13 which is considered as 38 strong turbulence. It is apparent from the plot that the weak theory greatly underestimates the aperture averaging values, and should not be used in considering such a strong turbulence strength. 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 5 10 15 Aperture Diameter (in cm) A p er t u r e A ver ag i n g F act o r , F Strong Theory Weak Theory C n 2 = 1e-13 (? I 2 = 4.32) Fig. 2.12: Comparison of the Aperture averaging Factor using Fante?s weak and strong turbulence theory for a choice of variance = 4.32 which is considered as strong turbulence level. Instead Equation (2.33) should be used for C n 2 in the weak turbulence range and Equation (2.42) for C n 2 in strong turbulence. When such a comparison is done in Figure 2.13, it is clear that the stronger the turbulence, the sharper the initial drop causing a faster saturation or a longer tail. Figure 2.13 provides an accurate measure of the relationship between the Aperture Averaging Factor values in weak and strong turbulence regimes. 39 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0 5 10 15 Aperture Diameter (in cm) A p er t u r e A v er ag i n g F act o r , F 1e-16 (4.32e-3) - Weak Theory 1e-13 (4.32) - Strong Theory C n 2 (? I 2 ) Fig. 2.13: Comparison of the Aperture Averaging Factor using Fante?s weak correlation functions with a variance of 4.32e-3 implying weak turbulence and a second curve using Fante?s strong correlation functions with a variance of 4.32 implying strong turbulence. 2.8 Conclusions This chapter presented an overview of the key turbulence theory that will be used throughout this dissertation. It also includes efficient computational techniques for Fante?s correlation functions that are important in assessing the effects of turbulence in weak and strong conditions. The aperture averaging factor was defined and presented for weak and strong turbulence using Andrews and Churnside approximations as well as the derived Fante?s correlation functions. 40 Chapter 3: Experimental Setup and Results 3.1 Introduction The foundations of the study of atmospheric turbulence were laid in the late 1960s and 1970s. Even after several decades of study, inconsistencies remain in the application of atmospheric turbulence theories to experimental systems, and the demonstration of acceptable agreement with experimental results. Theory reliably describes the behavior in the weak turbulence regime, but theoretical descriptions in the intermediate and strong turbulence regimes are less well developed. Most experiments to 1970 have been conducted over very long paths and were susceptible to the effects of scintillation saturation [28]. Other experiments done over shorter paths to exclude the saturation effects were inaccurate and inefficient in calculating the scintillation effects. By 1991, Churnside had conducted an experiment over 100, 250, 500, and 1000m paths to avoid saturation effects [3, 6]. He produced results in reasonably good agreement with spherical wave theory over short paths of 250 m. Over longer paths though, the results diverged from theory. Shortcomings of the experiment include the lack of the measurement of the background light which must be subtracted from the data for accurate measurements. Also, the scintillometer used in the experiment measured turbulence over a 250 m path, and not the path length under test causing inaccuracy in the calculation of C n 2 . In addition, only six apertures were 41 used in the experiment ranging from 1mm to 5cm, where larger aperture sizes are used by many commercial free space optical communication (FSO) units, so more measurements are needed. Linda Wasiczko, formerly of the Maryland Optics Group at the University of Maryland, conducted an experiment to measure aperture averaging over the same 863 meters link that was used in this thesis [29]. Her earlier work differed from the work described here because of differences in the receiver and the method of capturing and analyzing data. Dr. Wasiczko used two receiver apertures: a point receiver and a variable aperture receiver. The point receiver (5mm in diameter) is a scintillometer and is used to calculate the path-averaged C n 2 measurements. The variable receiver is a 20cm planoconvex lens with aperture stops ranging from 1cm to 16 cm. The beam is then received by a photodetector and the signal is recorded and processed in LabVIEW [29]. The main shortcoming of this experiment is the inaccuracy in data collection and the length of time required to collect all of the required data. Data collection for each aperture is typically done in 5 to 15 minute intervals, and to acquire proper data over all aperture sizes and various turbulence levels, data was mostly collected over several days. This provides a lot of inaccuracy in the measurement as well as inefficiency in the measurement time. For this reason, new empirical approaches that provide accurate and efficient experimental measurements are needed to help in the development of new theories. One such approach will be discussed in detail in this chapter. 42 3.2 Aperture Averaging Experimental Setup and Methodology A flexible empirical approach will be demonstrated for improving link performance through image analysis of intensity scintillation patterns coupled with frame aperture averaging on a free space optical (FSO) communication link. Aperture averaging calculations are invaluable in receiver design. A receiver must be large enough to collect sufficient power and reduce scintillation effects at a given range, but must also be of practical size. An imaging system has been constructed for measuring the effects of atmospheric turbulence and obscuration on FSO links. A He- Ne laser beam propagates over a range of 863 meters in atmospheric turbulence conditions that vary diurnally and seasonally from weak to strong. A high performance digital camera with a frame-grabbing computer interface is used to capture received laser intensity distributions at rates up to 30 frames per second and various short shutter speeds, down to 62.5?s per frame. The captured image frames are analyzed in LabVIEW to evaluate the turbulence parameter C n 2 , temporal and spatial intensity variances, and aperture averaging. 3.2.1 Experimental Setup An experimental system has been constructed for measuring the effects of atmospheric turbulence and obscuration on line-of-sight laser communication links. A 21mW He-Ne laser propagates between 2 rooftops over a range of 863 meters of free space at an average height above the ground of 12 meters. The propagating beam has a beam divergence on the order of 1.15 mrad. Figure 3.1 shows the University of Maryland campus map which locates the link range used with the transmitter location 43 shown on the roof of the A.V. Williams Building and the receiver at the Chesapeake Building. Fig. 3.1: University of Maryland Map showing the 863 meters link range used with the Transmitter shown on the roof of the A.V. Williams Building and the receiver at the Chesapeake Building. The transmitter uses a 21mW JDS Uniphase HeNe laser operating in a single mode (TEM 00 ) at 632.8 nm, with a 0.70 mm output diameter and 1.15 mrad beam divergence. The laser is followed by two reflective mirrors that direct the laser propagation direction towards the receiver end. The laser beam is then passed through a 30x Melles Griot beam expander which is adjusted to give a beam diameter at the receiver of approximately 1 meter. The transmitter is placed on the roof of the A.V. Williams building which is approximately 14 meters above the ground. A photograph of the transmitter is shown in Figure 3.2. The laser light then propagates over an 863 meter free-space path between the A.V. Williams Building and the Chesapeake Building mostly over asphalt parking lots with a few trees and grassy fields. Tx Rx 44 Rx HeNe Laser Beam Expander Mirrors Fig. 3.2: Transmitter of the 863m link from AVW to Chesapeake Building. The receiver consists first of a Meade LX200 Schmidt-Cassegrain telescope. Light is then filtered through a laser line filter at 632.8 nm with a 10 nm passband width. The filter removes any stray light from interfering with the operation of the receiver. The laser light then goes into a Pulnix TM1400 digital monochrome camera to be processed. The receiver in the Chesapeake building is approximately 12 meters above the ground. A photograph of the receiver is shown in Figure 3.3, and a complete schematic of the aperture averaging setup is shown in Figure 3.4. The receiver collects part of the wavefront that has been transmitted down range through the turbulent atmosphere. If the receiver has a small collection area, then the variance of the intensity that it will see is determined by the range length L, and the turbulence level. If the area of the receiver is increased, then the intensity variance decreases. This is to be expected, as in the limit a sufficiently large detector will collect all the transmitted light, and no atmospheric-turbulence-induced intensity variations should be seen. 45 Schmidt-Cassegrain Telescope Tx CCD Camera Fig. 3.3: Receiver of the 863m link from AVW to Chesapeake Building. Incident He-Ne radiation To LabVIEW Schmidt-Cassegrain Telescope Camera Filter He-Ne Laser Incident He-Ne radiation 30x Beam Expander Transmitter Receiver He-Ne Mirror He-Ne Mirror Fig. 3.4: Aperture Averaging Setup using a Schmidt-Cassegrain telescope with a 406.4 mm outer aperture size and a 127 mm inner obstruction diameter. The incoming beam intensity distribution is directed to a CCD camera to measure the variance of the received irradiance. 3.2.2 Frame Analysis Aperture averaging at the receiver is performed through the use of a digital monochrome 1392x1040 pixel camera with 4.65x4.65 ?m 2 pixels capturing frames at a rate of 30 frames per second and various short shutter speeds. A Schmidt- Cassegrain telescope is used to reproduce the intensity fluctuation pattern on the front of the telescope onto the camera CCD array, in a manner similar to that described by Moore, et. al. [30]. The telescope has an effective focal length of 4 meters at the 46 output flange of the telescope with a 406.4 mm outer aperture size and 127 mm inner obscuration diameter resulting in a total annular collection area of 1170 cm 2 . Additional imaging elements are not required since the telescope is setup in normal adjustment where the object is a collimated beam resulting in a demagnified output beam to the CCD. The captured frames are then analyzed in LabVIEW to calculate the turbulence index C n 2 , temporal and spatial variances, and the aperture averaging factor. Figure 3.5 shows a typical captured frame, with the characteristic multiple ?blob?-like variations of intensity across the field of view that is characteristic of laser beams significantly affected by turbulence. The central obscured region of the Cassegrain system is clearly visible. This image is essentially a compressed version of the input family of almost parallel rays entering the telescope front aperture. Fig. 3.5: Captured frame using a shutter speed of 1/500 second showing turbulence scintillation effects. The following steps describe the frame analysis performed on the image. First, the CCD camera settings are adjusted to acquire N number of images at a particular 47 shutter speed and frame rate. An AVI file is then acquired through a LabVIEW code interfacing with the camera. The file will contain a series of image frames which contain all of the frames information including the frame rate and the number of captured frames. Each captured frame is then subtracted from the background light. The background light is calculated through capturing a frame with the laser beam turned off. The resulting image frames can then analyzed for the required turbulence parameters. For a given range/wavelength scenario, the intensity variance for a point receiver is first measured and then the turbulence parameter C n 2 can be assessed through Equation (2.20). In practice, for a given range and wavelength, a receiver whose diameter is much less than the Fresnel zone size, D< () ()([].2/2/ 2 1 2 1 IiPIiP ppBER NN zeroone >+?<= += ) (4.13) Using Equation (4.11), the BER can be written as, . 22 1 2 1 22 1 2 1 22 1 2 1 2 1 22 1 2/ 2/ 2/2/ 2222 ? ? ? ? ? ? ? ? = ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + ? ? ? ? ? ? ? ? = ? ? ? ? ? ? += ?? ? ?? ? ?? N S erfc N S erfc N S erfc dxedxeBER I I xx ?? ?? (4.14) 71 In the presence of optical turbulence, the probability of error in Equation (4.14) must be averaged over the intensity fluctuation corresponding to a received ?one?. In this case, at any instant the actual received S/N ratio will rise or fall depending on whether the turbulence causes a ?fade? or a ?surge?. The probability of a ?zero? error does not change because there is no received signal when a ?zero? arrives and can be written as [32], ? ? ? ? ? ? ? ? = N S erfcp zero 22 1 2 1 . (4.15) However, the probability of a ?one? error is now, ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?= N S ierfcp one 2 1 2 1 2 1 , (4.16) where i is the detector signal corresponding to the actual received normalized intensity. To compute the average BER, the product of the probability of a ?one? error and the probability distribution for the normalized intensity (Equation (4.10)) must be integrated over all possible intensity values [32], () ( ? ? ? ? ? ? ? ? + ? ? ? ? ? ? ? ?= 0 2 2 22/1 2 )2/1(ln 2 1 exp 1 2 1 2 1 2 1 2 1 dii iN S ierfcp l l l one ? ? ?? , (4.17) which gives, ( zeroone ppBER += 2 1 ). (4.18) Figure 4.2 shows calculations of the average BER versus S/N ratio in the absence of turbulence. The threshold is set to 1/2. Different curves are plotted for 72 different values of intensity The plot shows the expected reduction in BER as the variance is decreased, implying better system performance at lower turbulence levels. -14 -12 -10 -8 -6 -4 -2 0 0102030 S/N Ratio (dB) lo g < B E R > 0 (No Turbulence) 0.23 0.29 0.327 0.375 0.467 0.736 1 Variance (? I 2 ) Fig. 4.2: Average BER versus S/N ratio in the absence of turbulence. Threshold is set to 1/2, and different curves are plotted for several intensity variances ( ). 2 I ? A desirable BER for an FSO link is 10 -9 , and Figure 4.2 allows the base receiver S/N ratio required to achieve this value to be assessed for different received signal intensity variances. For example, for a received intensity variance of 0.29, the base receiver S/N required is about 26dB. It has been shown that for very high levels 73 of turbulence, additional improvement in BER can be achieved by reducing the threshold decision level for ?ones? and ?zeros? below the usually accepted value of 1/2 [32]. 4.4 Effective Signal-to-noise Ratio analysis The discussion to follow will present the relationship between a so-called effective S/N ratio and the receiver aperture size. The effective S/N ratio will be defined in detail in this section. However, it is important to note from the start that such a parameter is not equivalent to the S/N ratio defined in the previous section except for the no turbulence case where they become equivalent. This means that the effective S/N ratio can not be used to predict the BER as shown for the S/N ratio in Figure 4.2. This is due to the fact that the effective S/N ratio fails to show the expected increase in value when the transmitted intensity is increased. In order to account for higher transmitted intensity, an adaptive method in finding the optimum threshold level must be included in the calculation. Nevertheless, such a parameter serves to guide the appropriate selection of the receiver size that ensures optimum performance. Such a value will be compared with the experimental choice of optimum receiver size presented in the previous chapter. Under weak turbulence conditions, the radial intensity distribution is the one appropriate to a Gaussian beam [4], , 2 exp),( 2 2 2 2 0 0 ? ? ? ? ? ? ? ? ?= EE w r w w ILrI (4.19) 74 where r is the lateral distance from the center of the beam, w o is the minimum spot size at the transmitter, ( )[ ] 2/1 6/5 22 /233.11 LILE kwLww ?+= is a measure of the effective turbulence induced beam spot size, 222 0 ?Lww L += is the spot size of the phase front at range L in the absence of turbulence, and ? is the beam divergence angle. The transmitted intensity is defined as, , where P 2 0 /2 wPI txo ?= tx is the transmitted Power. For our experimental setup defined in the previous section, L = 863 m, ? = 632.8nm, P tx = 20 mW, w o = 15 mm, and ? =?3 milliradians. In the presence of atmospheric turbulence between the transmitter and receiver, the average received signal power can be approximated by [8], () ?? ?= ? ? ? 2 0 2/ 0 2 ,,0 8 ),( D s LIDrdrdLrIP (4.20) where D is the diameter of the receiver aperture. It follows that the mean signal current is , ss Pi ?= where ? is the responsivity of the photodetector. The major sources of noise in FSO communication are shot noise and Johnson noise. Shot noise originates from fluctuations in the rate of photo-produced carriers, with mean zero and variance sSN ieB2 2 =? . Johnson noise results from thermal fluctuations in receiver electronics, also with mean zero and variance . Since these noise contributions are statistically independent, the total noise has zero mean with variance RkTB JN /4 2 =? . 222 JNSNN ??? += The total noise is assumed to be additive, which means the output current from the detector is a random variable Ns iii += , with mean s i and variance 75 . 222 NSi ??? += The averaged signal power to noise power ratio which we call the effective signal-to-noise ratio is given by [8], , .. 2 33.11 0 2 6/5 2 2 0 2 2 SNRF kw L SNR i SNR I L I i s eff ?? ? + ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? + == (4.21) where is the S/N ratio in the absence of the optical turbulence (receiver S/N ratio), and F is the aperture averaging factor. 22 0 / NS iSNR ?= Figure 4.3 plots the effective S/N ratio versus aperture diameter for several variances. In evaluating the Johnson and shot noise, a bandwidth B of 10 7 Hz, temperature T of 300K, and a resistance R of 50 ohms were used. The plot shows the improvement of S/N ratio with reduction in variance as expected. The experimental results for the aperture averaging factor for a variance of 1.38 described in the previous chapter were used in Equation (4.21) to plot the experimental curve shown. Figure 3.13 showed a higher aperture averaging factor for Fante?s strong theoretical curve than the experimental results which were both evaluated at the same variance of 1.38. Figure 4.4 is consistent with such analysis since the aperture averaging factor is shown in the denominator of Equation (4.21), so a higher value represents a lower effective S/N ratio. It is also clear from the experimental curve in Figure 4.3 that the S/N ratio saturates around an effective S/N ratio of 10 dB. Such an optimum value occurs at an aperture size of around 7cm which is the proved experimental optimum receiver size for the 1.38 variance data in Chapter 3. 76 0 10 20 30 40 50 0 5 10 15 Aperture Diameter (in cm) E f f ect i ve S i gn al - t o- N o i s e R a t i o ( i n d B ) 0 (No Turbulence) 0.0432 (Fante's Weak Theory) 0.432 (Fante's Weak Theory) 1.38 (Fante's Strong Theory) 1.38 (Experimental) Variance, ? I 2 Fig. 4.3: Effective S/N Ratio as a function of Aperture Diameter for an experimental variance of 1.38 plotted along with several theoretical variances calculated using Fante?s Correlation functions. We can redefine Equation (4.21) in terms of the aperture-averaged variance through the use of Equation (2.56). In addition, if we assume that turbulence-induced beam spreading is small, the second parameter in the denominator of Equation (4.21) can be removed leading to, 77 () , .1 , 0 2 02 0 SNR SNR SNRSNR AA AAeff ? ? + = (4.22) where is the aperture-averaged variance. 2 AA ? 0 1020304 Receiver Signal-to-Noise ratio (dB) 0 -10 0 10 20 30 40 Eff e ct ive Signal- to-Noi se r a ti o (dB) EFFECTIVE SIGNAL-TO-NOISE RATIO no turbulence AA variance 0.01 AA variance 0.05 AA variance 0.1 AA variance 0.5 AA variance 1 AA variance 1.5 for various aperture-averaged (AA) variances Fig. 4.4: Effective signal-to-noise ratio for various aperture-averaged variances. Figure 4.4 plots the effective S/N ratio for various aperture-averaged variances. In the absence of turbulence, the effective S/N ratio is equivalent to the receiver S/N ratio as expected. In the presence of turbulence, the higher the aperture- averaged variance the lower the effective S/N ratio. Equations (4.21) and (4.22) fail to 78 show the expected increase in the effective S/N ratio when the transmitted intensity is increased. In order to account for higher transmitted intensity, an adaptive method in finding the optimum threshold level must be included in the calculation. Figure 4.4 however serves as a good comparison of the results in Figure 4.3. For the experimental aperture-averaged variance of 0.138 at the optimum diameter of 7cm, Figure 4.3 shows an effective S/N ratio of around 10dB. Figure 4.4 also shows the maximum effective S/N ratio at an aperture-averaged variance of 0.1 to be approaching 10dB. Such plots serve to guide the appropriate selection of receiver size that ensures optimum performance for a particular turbulence level. From Figures 4.2, 4.3 and 4.4, we can see the reduction in BER and increase in the effective S/N ratio with the reduction in the variance. This shows the improvement in the overall system performance with decreasing variance. Therefore, aperture averaging can significantly improve the performance of the link, especially as the turbulence gets stronger because of the faster saturation of the effective S/N ratio, as shown in Figures 4.3 and 4.4, which leads to a smaller choice of optimum receiver size for best performance in the face of scintillation. 4.5 Conclusions Analysis of BER and S/N ratios for different turbulence levels shows that aperture averaging can significantly improve the performance of the link, especially as the turbulence gets stronger. The data presented is valuable in guiding the design of receivers for FSO communication systems. This is especially true for agile FSO transceivers, where size and weight compromises are needed. It is important to design such systems with apertures that are large enough for satisfactory 79 performance, but where excessively large receiver apertures, which provide only a marginal improvement in intensity scintillation, are avoided. 80 Chapter 5: Geometrical Simulation Models 5.1 Introduction Atmospheric turbulence has a significant impact on the quality of a laser beam propagating through the atmosphere over long distances. Turbulence causes the optical phasefront to become distorted from propagation through turbulent eddies of varying sizes and refractive index. Turbulence also results in intensity scintillation and beam wander, which can severely impair the operation of target designation and free space optical (FSO) communications systems. A new model is presented in this chapter to assess the effects of turbulence on laser beam propagation in such applications. The atmosphere is modeled along the laser beam propagation path as a spatial distribution of spherical bubbles. The size and refractive index discontinuity represented by each bubble are statistically distributed according to various models. For each statistical representation of the atmosphere, the path of a single ray, or a bundle of rays, is analyzed using geometrical optics. These Monte Carlo techniques allow us to assess beam wander, beam spread and phase shifts along the path, as well as aperture averaging effects at the receiver. An effective C n 2 can be determined by correlating beam wander behavior with the path length. This model has already proved capable of assessing beam wander, in particular the (Range) 3 dependence of mean-squared beam wander, and in estimating phase shifts developed across the laser phasefront as it propagates 81 through turbulence as well as aperture averaging at the receiver. The Monte Carlo simulations are compared and show good agreement with the predictions of wave theory [33]. A Random Interface Geometric approach which closely resembles this technique has already been developed by the Maryland Optics Group at the University of Maryland and is presented in detail in [34]. Such a Random Interface approach models the extended random medium as a series of random curved interfaces with random refractive index discontinuities across them [34].The new geometrical approach presented in this chapter which models the random medium as a spatial distribution of spherical bubbles is equivalent in analysis to the Random Interface Model in simulating beam wander and phase shifts. However, the bubble model allows for control of separation between the different bubbles, which better describes the real eddy structure formation of the atmosphere. In addition, the Random Interface Model propagated just a single ray or a pair of collimated pencil- thin rays through the random medium. In the bubble model, a bundle of rays that are initially Gaussian or uniformly distributed, are propagated through the simulated random medium. A varying aperture size is also added at the receiver to measure the aperture averaging effect. This is the first geometrical model to date to simulate the aperture averaging effect. 82 5.2 Theory Overview 5.2.1 Various Approaches The theory of wave propagation in random media is contained in a rich literature based on different approaches to the problem. These approaches range from rigorous diffraction theory to completely physically-based heuristic theories. Interestingly, heuristic theories often yield similar results to the diffraction theory, and give a much more intuitive feel to the problem. In all these methods, the polarization of the wave is not considered, since polarization changes of the wave as it propagates through random media have found to be negligible [34,35]. The earliest attempt to define meaningful field statistics for a random medium was based on geometrical optics. The best known published work that illustrates this approach is Chernov?s book [21]. Chernov based his analysis on Fermat?s principle in which each ray takes the shortest optical path. The problem then simplifies to a standard variational problem, that of solving the Lagrangian to find the trajectories of the rays. The main limitation of this approach is that it neglects diffraction effects [34]. To overcome this limitation and obtain a more realistic theory, the wave equation may be solved by the method of small perturbations [35], called the ?classical? theory of optical propagation because it was the first method to derive its results directly from Maxwell?s equations. The major advantage of this technique is the transform of the original Helmholtz equation, which is a homogeneous partial differential equation with random, space-dependent coefficients, into an equation which has constant coefficients and a source term [34]. 83 An alternate approach to solving the statistics of the perturbations on a wave traveling through a random medium involves using a spectral phase-screen technique, first proposed by Lee and Harp [36]. Contrary to the method of small perturbations, the technique is physical and geometrical rather than mathematical, and any approximations which are made are done in a physical context, making the implications of the approximations easier to assess. Furthermore, the technique takes diffraction effects into account with a minimum of mathematics, and the final results are identical to the method of small perturbations. The Random Interface Model closely resembles this technique and is presented in detail in [34]. 5.2.2 Phase Structure Function We can show that a heuristic derivation of the five-thirds phase structure function behavior can be used in order to illustrate how a simple model can yield results similar to more rigorous theory [35, 37]. Consider two parallel rays separated by a distance r, and traveling from a source at z = 0 to a receiver at z = L. Clearly, those eddies having the most impact on the phase difference between the two beams are on the order of r. Eddies much smaller than r will tend to be uncorrelated between the two beams and hence average to zero. Larger eddies will tend to contribute equally in phase to both beams and thus will not affect the phase difference. Thus, we may divide the path traveled by both rays into N segments of length r, where N = L/r. The phase difference )(r?? due to a single segment between the two rays is then simply, ),()( rnkrr ???? (5.1) 84 where is the difference in index of refraction between two points separated a distance r, and k is the wave vector. Since by assumption )(rn? )(rn? is zero, )(r?? is zero. The phase structure function for a single segment is given by, )()( 2222 rnrkr ???? . (5.2) But the last term of Equation (5.2) is simply the refractive index structure function, which for a turbulent region is described by Equation (2.5). Thus, the phase structure function due to a single element is, .033.0)( 38222 rkCr n ??? (5.3) Finally, if we assume that the phase difference over different segments is statistically uncorrelated, we can average over the entire path of N segments to obtain the phase structure function, .033.0)()( 35222 LrkCrNrD np ??? ? (5.4) Tatarski has derived the phase structure function for a plane wave to be [27], ( ) 3/522 64.0 LrkCrD np = . (5.5) If instead of a plane wave, we considered a pair of collimated parallel pencil thin beams separated by a distance r, Equation (5.5) would exactly result (except for a factor of 2, ( ) 3/522 32.0 LrkCrD np = . (5.6) 5.2.3 Beam Wander In weak turbulence a laser beam will retain its beam shape, but the beam will be additionally broadened with distance traveled by turbulence over and above its 85 natural broadening. In addition the beam centroid will wander about. If a photograph of the beam is taken, both these effects will broaden a time-exposed picture of the beam cross-section. The beam wander is generally more important than the turbulence induced beam broadening for many optical communication links. In such links when a wide, fairly well collimated laser beam is used, a good approximation to the mean-squared beam wander is defined by Ishimaru as [11], ()( )[] 33/1 0 2 2 2 2 1 2 02 2.21 2 zlCzz W nl ? +?+= ??? , (5.7) where z is the propagation distance, is the inner scale, 0 l 02 1 R=? , where is the radius of the equivalent Gaussian wave, is the spot-size, and 0 R 0 W ( ) 2 01 W??? = . For plane waves, the equation for the mean square beam wander simplifies to [35], 33/1 0 22 2.2 zlC nl ? =? . (5.8) At a range of 1km with mm1 0 =l and , the root-mean-squared (RMS) beam wander is about 5mm. With , which would strictly speaking not correspond to ?weak? turbulence, the RMS beam wander would be 50mm. 3/2152 m10 ?? = n C 132 10 ? = n C If the transmitted beam spot size is , then the beam diameter at range L is 0 w .22 0 0 w w L w += ? ? (5.9) For example, with mm20 0 =w , m3.1 ?? = , and L = 1km, . With the received beam diameter increases to 416mm. In this chapter, we will deal only with the plane wave case: a large diameter Gaussian beam that is diverging mm402 =w mm1 0 =w 86 slowly is equivalent to a plane wave provided w< 0 L m3.1 ? , for example, beam breakup effects might not occur until L > 10,000km. This prediction is clearly incorrect since beam breakup depends on turbulence strength and has been observed over path lengths on the order of 1km. Since the atmosphere fluctuates on a time scale that goes up to around 1kHz, a fast framing camera is required to visualize the actual breakup of the beam. Experiments performed with point detectors or with detectors that perform aperture averaging will not reveal beam breakup effects directly, but will fold this effect into the statistical variations observed. It is probably safe to assume that if weak turbulence theory does not hold for a given range, ? and , that beam breakup might occur. When beam breakup occurs the concept of beam wander becomes less relevant [35]. 2 n C 5.3 The Spherical/Bubble Model The Spherical Bubble Model assesses the effects of turbulence on laser beam propagation. The atmosphere along the laser beam propagation path is modeled as a spatial distribution of spherical bubbles. The size and refractive index discontinuity represented by each bubble are statistically distributed according to various models. For each statistical representation of the atmosphere, the path of a single ray, or a bundle of rays, is analyzed using geometrical optics. These Monte Carlo techniques 87 have already proved capable of assessing beam wander, in particular the (Range) 3 dependence of mean-squared beam wander, and in estimating phase shifts between rays as the laser beam propagates through turbulence. An effective C n 2 can also be determined by correlating beam wander behavior with the path length. In addition, this model is used to simulate the aperture averaging effect at the range for choices of varying receiver aperture sizes. This is the first geometrical simulation to predict the aperture averaging factor. 5.3.1 Beam Wander Simulation The beam wander simulation calculates the three dimensional trajectory for a single ray traveling a distance L through a simulated random medium. x z y Target Length L Laser Beam Snell?s Law (3 dimensions in simulation) n1 sin( ? 1) = n2 sin( ? 2) ?1 ?2 n1 n2 Fig. 5.1: Spherical/Bubble Model - Beam Wander Simulation. 88 The random medium is modeled as a series of random spheres/bubbles of random refractive index. For each statistical representation of the atmosphere, the path of a single ray is analyzed using geometrical optics. The mean square beam wander is averaged over each run. Figure 5.1 illustrates the beam wander simulation. 5.3.1.1 Simulation Procedure and Geometrical Analysis The uniform and Gaussian (normal) distributions will be used in the simulation model to evaluate the refractive indices of the spheres as well as the separation between the spheres. For this reason, their definition is included in this section. A uniform distribution in the range of [0, 1], which we call here U(0,1) simply returns a number within the specified range with equal probability of any of the numbers in between to be chosen. In order to modify such a distribution to increase the range, a multiplicative constant needs to be applied to the uniform distribution, U(0,1). In addition, the end points can simply be modified by adding and subtracting the appropriate values to the random uniform distribution. A Gaussian (normal) distribution is one that has the following density distribution, ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ?= 2 2 1 exp 2 1 )( ? ? ?? x xf for ????? x , (5.10) where ? is the mean, ? 2 is the variance, and X here is said to possess a normal distribution which we denote here as N(?,? 2 ). A normal distribution can be obtained using several methods, one of which is described in detail by Morgan [38]. Only the final result is presented here and will be used in the following simulations. Morgan 89 has defined a method of forming a normal distribution from two uniform distributions. If we call the uniform distributions, V 1 = U(-1,1) and V 2 = U(-1,1), then we can generate a pair of independent normal variables with zero mean and unit variance from the relation [38], 21 21 ln2 ? ? ? ? ? ? ? = W W VN , (5.1) 21 21 ln2 ? ? ? ? ? ? ? = W W VN , (5.12) where . To generate a normal variable with mean ? other than zero, we can simply add ? to N 2 2 2 1 VVW += 1 or N 2 . To create a normal variable with standard deviation ?, we need to multiply the standard normal by ?. A three-dimensional Snell?s law needs to be evaluated in order to calculate the refracted output vector given the incident vector, the normal to the sphere, and refractive indices of spheres. Such analysis is described in detail below since it will be used in all of the simulations described in this chapter. When a ray traveling in free-space encounters a dielectric interface, it refracts due to the change of index it encounters. The boundary conditions state that the tangential component of the electric field and the normal component of the magnetic field are continuous. In order to satisfy such conditions, Snell?s law defines the relation between the incident and refracted angle as follows, )sin()sin( 2211 ?? nn = , (5.13) where n 1 and n 2 are the indices of refraction of air and the sphere respectively when the ray is entering the sphere and reversely when the ray is exiting the sphere. ? 1 and 90 ? 2 are the incident and refracted angles with the normal to the spheres, respectively. In addition, the boundary conditions require that the refracted ray 2 r r to lie in the same plane as the plane of incidence which includes the incident ray and the normal to the sphere . We define 1 r r n r 1 P r to be the perpendicular vector to the plane of incidence with magnitude ( ) 1 sin ? and 2 P r to be in the same direction as but with magnitude . If we use unit vectors for the incident ray ( ), refracted ray ( ) as well as the sphere?s normal ( ), then, we can write 1 P r ( 2 sin ? ) 1 ?r 2 ?r n? 1 P r and 2 P r as, nrP ?? 11 ?= r , (5.14) nrP ?? 22 ?= r . (5.15) Given the incident ray and the sphere?s normal , 1 ?r n? 1 P r can be determined through Equation (5.14). The incident angle is then 1 ? , ( ) 11 sin P r =? . (5.16) The refracted angle 2 ? can then be calculated through the application of Snell?s law, Equation (5.13). Therefore, 2 P r is, ( ) 212 sin ? ??= PP r , (5.17) where 111 ? PPP rr = is the normalized unit vector of 1 P r . 91 Figure 5.2 show the three-dimensional Snell?s Law which can be used to determine the refracted ray . The refracted unit vector can be evaluated through adding the 2 ?r 2 ?r 2 ? Pn r ? vector which is in the plane of incidence to a scaled version (by a constant ? ) of the unit vector of the sphere?s normal , n? 22 P?n ? ?? r ?+= nr . (5.18) Therefore, in order to calculate the output vector , we only need to know the incident vector, the normal to the sphere, and the refractive indices of the spheres assuming air of index 1.0001 in between them. 2 ?r ? 2 ? 2 ? 1 n 1 n 2 n ? ? 1 r r 2 ? r? (wrong ray) n ? ? 2 r r 2 ? Pn r ? 2 ? r R Fig. 5.2: Three dimensional Snell?s Law. Equation (5.18) can be rewritten as, 2 2 2 2 2 2 ??2??? PnnnPnr rr ??++?= ?? , (5.19) 92 where 1? 2 =r .But the last term in Equation (5.19) is zero since the vector is perpendicular to the n? 2 ? Pn r ? vector as shown in Figure 5.2. If we specify [cban ,,P? 2 =? ] r and [ ] zyx nnnn ,,? = , then . (5.20) 01 2222 =?+++ cba? Defining , and , and noting that , the solution for ? can be written as, 1a ~ = 0222b ~ =++= zyx cnbnan 1c 222 ~ ?++= cba 1nnn 2 z 2 y 2 x =++ ~ ~~ ~ 2 ~ a2 ca4bb- ? ?? = . (5.21) This equation yields two solutions for , but only the solution that gives a positive z component is correct and is chosen. Using Figure 5.2, the unit vector can be multiplied by a scaling factor to determine the full refracted vector , 2 ?r 2 ?r 2 r r ( )( ) 222 cos2? ?Rrr ?= r , (5.22) where R is the radius of the sphere. A particular simulation run is comprised of the following: The user inputs the length of the target L, the mean free path l ? which defines the distance traveled in free-space before the sphere is encountered, and standard deviation for the path l ? , the mean index of refraction of the spheres n ? taken as the air index and standard deviation for the refractive index n ? , and the number of runs desired N. A beam is started out at location (0,0,0) in the z-direction. The index of refraction of the spheres 93 is chosen from a Gaussian distribution, and the sphere?s center and radius are chosen from a uniform distribution with the condition that the propagating ray must hit the chosen sphere. Snell?s law is then used to evaluate the new output vectors at the entering and exiting point of the sphere. The ray then travels through free-space a certain path length, chosen from a Gaussian distribution with mean l ? and standard deviation l ? . The rays then encounter another sphere. This entire process is repeated until the target length is reached. The mean square of the distance from the center of the target 2 l ? is updated, and the next run begins. In this simulation, no correlation between spheres is assumed through choosing the mean free path to be much larger than the inner scale of turbulence. BEGIN INPUT DATA WHILE RUN #< N END FALSE TRUE START RAY IN Z-DIRECTION AT (0,0,0) WHILE LTOTAL < L FALSE UPDATE TRUE PICK INDEX OF SPHERE FROM GAUSSIAN PICK SPHERE?S CENTER AND RADIUS FROM UNIFORM CALCULATE UNIT NORMAL OF SPHERE USE 3D SNELL?S LAW TO GET BEAM?S REFRACTED RAYS PICK LENGTH TO SPHERE l FROM GAUSSIAN PROPAGATE BEAM A DISTANCE l Fig. 5.3: Flowchart of the Beam Wander Simulation Model. 94 5.3.1.2 Beam Wander Simulation Results Figure 5.4 plots the beam wander simulation results using the following parameters: L = 1km, l ? = 1m, l ? = 0.9m, n ? = 1.0001 (free-space index), n ? = 0.00001, and N = 1000. The radius of the randomly selected spheres is chosen in the range of 1mm to 1m. It is clear that the simulation results show excellent agreement with the cubic fit described in theory in Equation (5.8). 0 200 400 600 800 1000 -0.002 0.000 0.002 0.004 0.006 0.008 0.010 0.012 0.014 M S B eam W ander (in m 2 ) Range (in m) Mean Square Beam Wander Beam Wander_Cubic Fit Fig. 5.4: Mean Square Beam Wander and Cubic Fit for a 1km range using the Spherical/Bubble Model. The nominal value for is evaluated through the use of Equation (5.8) and for this particular run was 5.87?10 2 n C -13 . This is a fairly strong turbulence level. Different values of can be obtained by varying the simulation parameters. For 2 n C 95 lower turbulence levels which can give a lower beam wander at the receiver, a variation in the simulation parameters can easily achieve this. For example, the choice of the bubbles? standard deviation of their refractive index fluctuations n ? can be reduced or the mean free path l ? can be increased to achieve lower turbulence levels and consequently decrease the beam wander. In addition, only the range can be reduced in the chosen simulation parameters for lower turbulence levels. Another plot is shown in Figure 5.5 where the range is chosen to be 5 km. The simulation parameters used were: L = 5km, l ? = 1m, l ? = 0.9m, n ? = 1.0001 (free- space index), n ? = 0.00001, and N = 1000. 0 1000 2000 3000 4000 5000 0.0 0.4 0.8 1.2 1.6 M S Be a m W a nder ( i n m 2 ) Range (in m) MS Beam Wander Cubic fit Fig. 5.5: Mean Square Beam Wander and Cubic Fit for a 5km range using the Spherical/Bubble Model. 96 Again, the mean square beam wander in Figure 5.5 follow the expected cubic fit described by theory in Equation (5.8). Since the range has been increased from 1km in the simulation plotted in Figure 5.3 to 5km in Figure 5.4, the turbulence level has also increased. This is shown from the increase in which is evaluated for this simulation to be 5.9?10 2 n C -13 . 5.3.2 Phase Wander Simulation Model The trajectories of two parallel rays propagating through the simulated random medium are computed simultaneously. Figure 5.6 illustrates the phase wander simulation procedure. The simulation is very similar to the beam wander simulation, except that two rays are now considered. Fig. 5.6: Spherical Bubble Model ? Phase Wander Simulation. The total distance traveled by each ray as both rays travel to the target is computed. The difference in the path length traveled l? will yield the phase difference ?? between the rays through [34], Beams Separated Distance r in mm x z Laser Beam y Target Length L Beam 1 Beam 2 l1 l1? l2 l3 l2? l3? l4? l5? l6? l7? l8? l9? l4 l5 l6 l7 l8 l9 l10 l11 l12 l13 l14 97 l?=? ? ? ? 2 . (5.23) The three-dimensional geometry of the sphere is taken into account in the propagation of the rays. The spheres are chosen such that the first ray intersects all of the selected spheres. Therefore, the first ray follows the same path described in the beam wander simulation model. The second ray though is started out in the positive z-direction separated a distance r from the first ray in the x-direction. First, a check is made to see whether or not the second ray intersects the selected sphere. This can be done by solving the sphere equation with the second ray vector equation. The line equation can be defined as, n zz m yy l xx 121212 ? = ? = ? , (5.24) where (x 1 ,y 1 ,z 1 ) and (x 2 ,y 2 ,z 2 ) are 2 points on the line, and l, m, n are the directional cosines of the line. The directional cosines are defined as l = cos?, m = cos?, n = cos?, where ?, ?, ? are the angles that the vector makes with the positive x-, y- and z- axes, respectively. In addition, we have the following inequality, . (5.25) 1coscoscos 222 =++ ??? The sphere equation is defined as, ()( ) ( ) 2 2 02 2 02 2 02 Rzzyyxx =?+?+? , (5.26) where (x 2 ,y 2 ,z 2 ) is a point on the sphere, (x 0 ,y 0 ,z 0 ) is the center of the sphere and R is the Radius of the sphere. A solution is checked for z 2 by plugging in x 2 and y 2 in terms of z 2 in Equation (5.24) into Equation (5.26). If a solution exists, the smaller root of z 2 is taken from which x 2 and y 2 can be solved for. Then, the three- 98 dimensional Snell?s law can be applied as in the procedure outlined in the previous section to calculate the refracted vector. If however no solution exists, this implies that the second ray can not intersect the selected sphere. In this case, the ray is made to propagate to the same z-location as the first ray. This can be achieved through solving the second ray?s incident vector with the plane z = k, where k is the z-location of the first ray after it propagated through the sphere. For example, let us assume the incident unit vector for the second ray to be (V x , V y , V z ), and its starting point as (x 1 , y 1 , z 1 ). Then the end point for the second ray is, ( )ktVytVxzyx yx ,.,.),,( 11222 ++= , (5.27) where . z Vzkt /)1( ?= A single simulation run is comprised of the following steps: The user inputs the length of the target L, the mean free path l ? , and standard deviation for the path l ? , the mean index of refraction n ? and standard deviation for the refractive index n ? , the starting separation between the beams r and the step size between separations ? r, and the number of runs desired N. Two beams are started out separated by a distance r at the coordinates (0,0,0) and (-r,0,0). The index of refraction of the spheres/bubbles for both beams is chosen from a Gaussian distribution. Similarly, the sphere?s radius and center are chosen from a uniform distribution with the condition that the first propagating ray starting at (0,0,0) must hit the randomly chosen sphere. A single path length is chosen from a Gaussian distribution, and the first beam travels this length to intersect the sphere. If the second beam can intersect the sphere, then it is made to propagate till it intersects it. Otherwise, if no intersection occurs, the 99 second beam is made to travel a path length such that the z-coordinates for both beams are made the same. Snell?s law is then invoked to find each beam?s new output vector. A running sum is kept of the distances traveled by each beam, and the entire process is repeated until the z-coordinate of each ray equals the target length. The square difference between the beams is then averaged and saved for the particular value of separation r, after which r is decremented by ? r and the simulation is restarted until r = 0 (for which 0=?? ). BEGIN INPUT DATA WHILE RUN #< N END FALSE TRUE START BOTH RAYS IN Z-DIRECTION AT (0,0,0) WHILE LTOTAL < L FALSE UPDATE TRUE PICK INDEX OF SPHERE FROM GAUSSIAN PICK SPHERE?S CENTER AND RADIUS FROM UNIFORM CALCULATE UNIT NORMAL OF SPHERE FOR BEAM1 USE 3D SNELL?S LAW TO GET BEAM1 REFRACTED RAYS PICK LENGTH TO SPHERE l FROM GAUSSIAN PROPAGATE BEAM1 A DISTANCE l PROPAGATE BEAM2 TO SAME Z AS BEAM1 WHILE SEPARATAION > 0 TRUE IF BEAM2 INTERSECTS SPHERE NO YES USE 3D SNELL?S LAW TO GET BEAM2 REFRACTED RAYS ? Fig. 5.7: Flowchart of the Phase Wander Simulation Model. 100 5.3.2.1 Phase Wander Simulation Results The mean square phase difference is taken to be the phase structure function, which was shown for a pair of collimated parallel pencil thin rays in Equation (5.6) and obeys a five-third law assuming weak turbulence. In Figure 5.8, the simulation results are plotted for varying mean free path l ? on a log-log graph. 0.70.80.9123456789102030 10 11 10 12 10 13 10 14 < ?? 2 (r )> Separation r (cm) ? l ,? l 1,0.9 0.7,0.7 0.5,0.5 0.3,0.3 0.1,0.09 0.01,0.01 5/3rd slope Figure 5.8: Phase Wander Simulation Structure Function for Various Mean Free Paths. The fixed parameters assumed in Figure 5.8 are: target length L = 1 km, mean index of refraction 00001.1= n ? with standard deviation 000001.0= n ? , and starting beam separation r = 0.25m with ? r = 0.01. The beam separation was chosen to range 101 from cmL 2510 =? to 0 (no separation) assuming a He-Ne laser beam propagating with ? = 632.8nm. The simulation was averaged for 1000 runs. Figure 5.8 show that as l ? is decreased, the y-intercept, which represents a number proportional to increases as expected. This is because a decrease in , 2 n C l ? causes more collisions to occur and thus larger path differences between the beams which represents a larger turbulence level. Figure 5.8 agrees well with this prediction. Furthermore, it is found that the slope of the curve is increasing with higher l ? values or lower turbulence levels. For the highest value of m l 1=? plotted, the slope is approximately 4/3rd. It is expected for a high enough choice of l ? that causes weak turbulence, the slope should follow the five-third power level predicted by theory. Although the form of the refractive index fluctuations was assumed Gaussian rather than following a Kolmogorov spectrum, perhaps the random encounters with many spheres tends to blur the difference between the assumed functional forms of the refractive index fluctuations, leading to the same statistical result in the wave parameter. 5.3.3 Aperture Averaging Simulation Model This model calculates the aperture averaging factor for a number of rays, Gaussian as well as uniformly distributed, propagating through turbulence simulated spherical bubbles into a circular receiver of varying aperture size. The three-dimensional trajectory of each ray is analyzed using geometrical optics. The number of rays that reach the target length L, within the selected receiver aperture size are summed. The variance of the total rays within each aperture size is 102 then calculated over N simulation runs. Such variances are then normalized by the variance of the smallest chosen aperture size to evaluate the aperture averaging factor, . )0( )( 2 2 = = D D F I I ? ? (5.28) Fig. 5.9: Aperture Averaging Simulation using the Spherical Bubble Model. 5.3.3.1 Simulation Procedure and Geometrical Analysis The three-dimensional space is first filled with spherical bubbles of varying sizes ranging between 1mm and 1m. There is a coverage percentage chosen for the percentage of bubbles filling the three-dimensional space. The bubbles are chosen such that they do not intersect or touch one another. The starting x and y coordinates of the rays (x 0 ,y 0 ) are each selected from a random Gaussian distribution with means 0= x ? , 0= y ? , respectively and variance , where w is the beam spot 22 w=? Selected Receiver Aperture Sizes x y (0,0,0) z Ydim Range, L Bundle of ys ra -Gaussian distributed Xdim 103 size at the transmitter (starting point) taken as 20mm. The starting z-coordinate of all of the rays is at z=0, and the beam divergence angle, ? ? at the transmitter is chosen as 1 milliradians. 5.3.3.1.1 Three-Dimensional Gaussian Ray Modeling The starting three-dimensional unit vector for each ray needs to be determined before the rays can be propagated through the randomly simulated medium, i iii i i i S kzjyix S r r ? ?? ++ == ? r , (5.29) where 222 iiii zyxS ++= is the magnitude of the vector i r r , and i the index of the ray. z x y x 0 i y 0 i R i ? x i ? y i O O 1 Fig. 5.10: Position Figure showing the origins of the rays in the x-y plane. 104 Figure 5.10 is a position figure showing the origin ( )0,, 00 1 ii yxO = of the rays in the x-y plane, where 2 0 2 0 iii yxR += , (5.30) is the distance from the Origin ( )0,0,0=O coordinate to the starting point of the ray, O 1 . Hence, ,cos 0 i x R x i i =? (5.31) i y R y i i 0 cos =? . (5.32) Figure 5.11.a shows the angles of the Gaussian randomly distributed rays in the x-z plane where, ,tan 1 1 i i i x z x z =? (5.3) and from our choice of the initial divergence angle at the transmitter, we can say that, 0 0 tan ?? ?= w x i ix z . (5.34) Since all of the parameters in Equation (5.34) are known, we can determine i x z ? through which and can be determined. i x 1 i z 1 Similarly, Figure 5.11.b shows the angles of the Gaussian randomly distributed rays in the y-z plane where, 105 i ii i y z y w y z 2 2 0 0 tan =?= ?? . (5.35) Equation (5.35) can be used to determine and . i y 2 i z 2 x z z 1 i x 1 i x 0 i ? z x i i r 1 r z z 2 i i r 2 r y ? z y i y 2 i y 0 i Fig. 5.11a: Rays in the x-z plane. Fig. 5.11.b: Rays in the y-z plane Figure 5.12 shows the random three-dimensional Gaussian ray propagation vector i r r originating from the x-y plane and propagating in the x-y-z direction, where, i xi zi xx ?sin 1 == , (5.36) i yi zi yy ?sin 2 == , (5.37) i y i xii zzi zzz ?? sinsin 21 +=+= . (5.38) Therefore, kzzjyixr iiii i ? )( ?? 2121 +++= r . (5.39) 106 x z y i r 1 r i r 2 r i r r x 1 i z 1 i z 2 i y 2 i O 1 Fig. 5.12: Random Three-dimensional Gaussian Ray propagation vector i r r . 5.3.3.1.2 Simulation Procedure A typical simulation run is comprised of the following: First the user is prompted to input the X and Y dimensions as well as the range L of the three- dimensional space. The user then inputs the coverage percentage (Spheres Cover) of the spheres as well as the mean n ? and variance n ? of the spheres? index of refraction, which will be calculated from a random Gaussian distribution. Finally the user inputs the beam?s starting divergence angle 0 ? and beam waist w as well as the 107 number of simulation runs N. The space is first filled with bubbles according to the percentage of coverage chosen. The bubbles are then sorted by the z-coordinates of their centers. The rays all start at z = 0 and their x and y coordinates are selected from a random Gaussian distribution with mean 0= x ? and 0= y ? , respectively and variance as described above. The starting angle of each ray with the axes is determined using the procedure outlined above from which each ray?s starting unit vector can be defined. Each ray is checked for intersection with any of the bubbles and the bubble with the smallest z-coordinate is chosen for the ray to refract through it according to Snell?s law. Such a check is done until no more spheres can be intersected by the ray within the chosen three-dimensional space. Then, the ray is made to propagate to the target length L. The ray?s beam wander from the center 22 w=? 22 yx + is checked to whether or not it lies within some selected circular aperture sizes at the receiver. If so, it gets counted. The same procedure is repeated for each of the rays. Then the whole simulation is repeated N times for different rays? distributions. The variance of the total rays within each aperture size is then calculated over the simulation runs. Such variances are then normalized by the variance of the smallest chosen aperture size to evaluate the aperture averaging factor, F. 5.3.3.2 Simulation Results 5.3.3.2.1 Turbulence Strength First, one ray starting at the origin is propagated in the z-direction through the simulated randomly chosen spheres. Such a simulation is performed in order to 108 determine the level of turbulence for varying standard deviation of the spheres? index of refraction n ? . The chosen simulation parameters are: L=1km, Spheres Cover = 20%, n ? = 1.0001 and N = 1000 runs. The mean-squared (MS) beam wander (Equation (5.8)) is calculated for varying n ? values. C n 2 can then be calculated through the Rytov variance equation where the inner scale of turbulence l o which represents the smallest sphere size is taken as 1mm. Table 5.1 shows C n 2 for varying n ? values. It is clear from the C n 2 values presented that a n ? of 1e-5 represents a relatively strong turbulence level, while a n ? of 1e-7 already gives a relatively weak turbulence level. Such values of n ? representing weak and strong turbulence levels will be used as a limiting choice in the simulations to follow. ? n MS Beam Wander (m) C n 2 1.00E-04 2.35E+03 1.0669E-07 5.00E-05 1.85E+01 8.4114E-10 1.00E-05 3.47E-03 1.576E-13 1.00E-06 2.39E-05 1.0852E-15 5.00E-07 5.08E-06 2.3086E-16 1.00E-07 2.91E-07 1.3218E-17 1.00E-08 9.38E-08 4.2641E-18 Table 5.1: C n 2 for varying n ? values. A n ? of 1e-7 already gives a relatively weak turbulence level while a n ? of 1e-5 represents relatively strong turbulence. 109 5.3.3.2.2 Ray Number Reaching Circular Receiver Apertures The simulation was then performed for the Gaussian distributed rays using the procedure outline in the previous subsection for 1000 Rays. The parameters chosen in the simulation are as follows: L=1km, Spheres Cover = 20%, n ? = 1.0001, 0 ? = 1mrad, w = 20mm, and N = 1 run. First the simulation number, N was set to 1 and 1000 rays that are Gaussian distributed were propagated through the randomly simulated bubbles. Figure 5.13 shows curves for the number of rays that fit within circular aperture receivers of varying sizes. Different curves were plotted for bubbles with varying n ? along with the no bubbles case. On average, the space was filled with about 14,000 bubbles with the rays intersecting about 300 bubbles on average along its path. The results show the expected reduction in received number of rays with increasing n ? or increasing turbulence level as expected. However, such a conclusion is only expected for severely strong turbulence levels. Table 5.1 show that a n ? of 1e-5 already represents a relatively strong turbulence level yet the ray number in Figure 5.13 is very close in value to the no turbulence level then saturating at a bit lower value. This is due to the fact that higher turbulence levels do not necessarily imply that the rays will miss the targeted aperture radius especially for small aperture radii. It only means a higher beam wander which could be accounted for through outer rays bending inwards towards the center and being collected by the aperture. However, severely strong turbulence bent the rays sufficiently to mostly miss the selected aperture radii. 110 0 200 400 600 800 1000 0.00 0.20 0.40 0.60 0.80 Aperture Radius (m) R a y N u m b er No Bubbles 1.E-05 5E-05 1.E-04 ? n Fig. 5.13: Ray Number verses Receiver Aperture Radius for 1 simulation run of Gaussian Distributed rays. Various curves are plotted for bubbles of varying standard deviation ?std?, n ? along with the no bubbles case. The same simulation procedure was repeated for plane rays instead of Gaussian rays. The rays are started uniformly distributed in the x-y plane and propagating in the z-direction. The rays are chosen to initially fill a circular region in the x-y plane equivalent in diameter to the largest aperture size chosen at the receiver side. Figure 5.14 shows the ray number for 1000 rays over 1 simulation run. 111 0 200 400 600 800 1000 0.00 0.10 0.20 0.30 0.40 0.50 Aperture Radius (m) Ra y Nu m b e r No Bubbles 1.E-05 5E-05 1.E-04 ? n Fig. 5.14: Ray Number verses Receiver Aperture Radius for 1 simulation run of Uniformly Distributed rays. Various curves are plotted for bubbles of varying standard deviation ?std?, n ? along with the no bubbles case. It is clear from Figure 5.14 that a n ? of 1e-5 which represents a relatively strong turbulence level can give even a higher ray number than the no turbulence case for the lower aperture diameters. This confirms the analysis described for the Gaussian Ray model in Figure 5.14 which implied that stronger turbulence can lead to an equivalent or even higher ray number reaching the smaller aperture sizes but only up to a certain degree of turbulence level. 112 The Gaussian ray distribution failed to show the higher ray number for stronger turbulence levels at small aperture diameters as was shown for uniformly distributed rays. This could be due to the fact that the Gaussian distributed rays are initially diverging outwards which increases the possibility of ray miss at the receiver especially for small aperture sizes. 5.3.3.2.3 Aperture Averaging Factor ?F? Now, the aperture averaging factor ?F? can be calculated through using Equation (5.28). The sample variance form is used in calculating the variance of a particular aperture diameter over 100 simulation runs. It is then normalized by the mean squared to calculate the normalized aperture-averaged variance, 2 2 2 1 )( )( i j iij iI N x D ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? = ? , (5.40) where i is the aperture diameter index and j is the simulation run index. Equation (5.40) is then divided by the normalized variance of the minimum aperture diameter to evaluate the aperture averaging factor F (Equation (5.28)). Through plugging Equation (5.40) in Equation (5.28), the aperture averaging factor is calculated and plotted in Figure 5.15 for the Gaussian distributed rays with 0 ? = 1mrad and n ? =1e-7, 1e-6 and 1e-5 corresponding to weak, intermediate and strong turbulence levels respectively. 113 0.00 0.20 0.40 0.60 0.80 1.00 0.00 0.20 0.40 0.60 0.80 1.00 Aperture Radius (cm) F Fa c t or 1.00E-07 1.00E-06 1.00E-05 ? n Fig. 5.15: Aperture averaging factor (F) versus the aperture radius for Gaussian distributed rays. Different curves are plotted for n ? =0.00001, 0.00005, and 0.0001. Figure 5.15 shows that the higher the standard deviation of the bubbles? indices n ? , the higher the aperture averaging factor. This is expected since the ray number in Figure 5.13 for even relatively strong turbulence such as 1e-5 is lower than the no turbulence case even at small aperture sizes. This implies that the Gaussian ray model used gave the expected reduction in ray number reaching the receiver aperture with increasing turbulence level at all receiver aperture sizes. Hence the mean squared value of the number of rays in the denominator of the aperture-averaged variance in Equation (5.40) will be less for higher turbulence levels. In addition, the variance in the numerator of Equation (5.40) is expected to increase with increasing 114 turbulence level. These two effects both imply an increase in the F Factor with higher n ? values. Figure 5.16 plots the aperture averaging factor F versus aperture radius for uniformly distributed rays. Different curves are plotted for several n ? values. 0.00 0.20 0.40 0.60 0.80 1.00 0.00 0.40 0.80 1.20 1.60 2.00 Aperture Radius (cm) A p er t u r e A ver ag i n g F act o r , F 0.00001 0.00004 0.00008 ? n Fig. 5.16: Aperture Averaging Factor versus Aperture Radius for uniformly distributed rays. Different curves are plotted for several n ? values all within the strong turbulence regime. Figure 5.16 show that the higher the n ? value, the higher the aperture averaging factor. This is because all of the n ? choices in the plot present strong turbulence levels. This implies that the higher the n ? value, the lower the mean square value of the number of rays reaching the particular aperture diameter. Since 115 the mean square shows in the denominator of the aperture-averaged variance (Equation (5.40)), then this corresponds to a higher F factor. Figure 5.17 shows the F factor for uniformly distributed rays plotted for n ? = 1e-7, 1e-6 and 1e-5 corresponding to weak, intermediate and strong turbulence strengths correspondingly. It is clear from the plot that the higher the n ? value which implies higher turbulence level, the sharper the initial decline followed by a longer tail. This exact result is expected and was proved in the theoretical analysis as well as from the experimental results shown in Chapter 2 and 3. 0 0.2 0.4 0.6 0.8 1 0.00 0.50 1.00 1.50 Aperture Radius (cm) F Fa ct or 1.00E-07 1.00E-06 1.00E-05 ? n Fig. 5.17: Aperture Averaging Factor versus Aperture Radius for uniformly distributed rays. Different curves are plotted for several n ? values corresponding to different turbulence strengths. 116 This result is expected from Figure 5.17 using such a geometric ray model since a n ? of 1e-5 corresponding to relatively strong turbulence has given a higher ray number at small aperture diameters when compared to the no bubbles case as observed in Figure 5.14. If we think of weak turbulence such as a n ? of 1e-7 to follow a similar behavior to the no turbulence case at low aperture diameters, then the mean square value of the ray number for small apertures will be higher for the n ? of 1e-5. Such analysis confirms that the F factor for n ? of 1e-5 to be less at small aperture sizes when compared to a n ? of 1e-7. However, at larger aperture diameters, the fluctuation from the mean term in the numeration of Equation (5.40) becomes more noticeable for stronger turbulence levels dominating the increase in the mean square ray number in the denominator. This effect lead to a higher F factor for stronger turbulence levels which explains the long tail at larger diameters. 5.4 Conclusions A new geometrical model was presented to assess the effects of turbulence on laser beam propagation in the atmosphere. The atmosphere was modeled along the laser beam propagation path as a spatial distribution of spherical bubbles with refractive index discontinuities that are statistically distributed according to various models. For each statistical representation of the atmosphere, the path of a single ray, or a bundle of rays, was analyzed using geometrical optics. These Monte Carlo simulations have proved capable of assessing beam wander, in particular the (Range) 3 dependence of mean-squared beam wander, and in estimating phase shifts between the rays as they propagate through turbulence and in calculating the aperture 117 averaging effect. An effective C n 2 was also determined by correlating beam wander behavior with the path length. The Monte Carlo simulations were all compared and show good agreement with the predictions of wave theory. This is the first report to date to simulate the aperture averaging factor using geometrical ray analysis. The aperture averaging factor was modeled for Gaussian as well as uniformly distributed rays. The results assessed the reduction in scintillation with increasing aperture diameter. In addition, for the uniformly distributed rays, the aperture averaging factor F showed the sharper initial decline followed by a longer tail for higher turbulence strengths as predicted by the theoretical as well as experimental results. 118 Chapter 6: Conclusions and Future Research A flexible empirical approach was demonstrated for improving free space optical (FSO) communication link performance through image analysis of intensity scintillation patterns coupled with frame aperture averaging on FSO communication link. The aperture averaging results shown in Chapter 3 demonstrated the expected reduction in intensity fluctuations with increasing aperture diameter, and show quantitatively the differences in behavior between various strengths of turbulence. The reduction in scintillation with aperture size guides the selection of the optimum receiver aperture. Efficient computational techniques for Fante?s correlation functions that are important in assessing the effects of turbulence in weak and strong conditions were developed in Chapter 2 to compare with the experimental results. The experimental results presented in the intermediate turbulence region fitted between the weak and strong turbulence theory and show ?excellent agreement? with the expected behavior. Such results show significant improvement upon existing empirical data in accuracy and should be very valuable in the development of new theories in the intermediate turbulence regime. Theory reliably describes the behavior in the weak turbulence regime, but theoretical descriptions in the intermediate and strong turbulence regimes are less well developed. For this reason, such results should help in the development of new theories that can help in filling the existing gap. 119 Spatial and temporal variance analyses within single frames and between frames were also compared and show good agreement. It was also observed that turbulence is weaker during light to moderate rain conditions due to the reduction in the intensity scintillation. This causes a better performance of the FSO links during rain which is one of its advantages over RF technology. In addition, during sunrise and sunset where the ground surface and surrounding air temperatures are at an equilibrium state, the lowest turbulence level within the day is noticed. It was finally proved that the irradiance calculations using the frame analysis are independent of the shape of the receiver aperture but instead only on the area. Analysis of BER and S/N ratios for different turbulence levels shows that aperture averaging can significantly improve the performance of the link, especially as the turbulence gets stronger. The data presented is valuable in guiding the design of receivers for FSO communication systems. This is especially true for agile FSO transceivers, where size and weight compromises are needed. It is important to design such systems with apertures that are large enough for satisfactory performance, but where excessively large receiver apertures, which provide only a marginal improvement in intensity scintillation, are avoided. A new geometrical model was also presented in Chapter 5 to assess the effects of turbulence on laser beam propagation in the atmosphere. The atmosphere was modeled along the laser beam propagation path as a spatial distribution of spherical bubbles with refractive index discontinuities that are statistically distributed according to various models. For each statistical representation of the atmosphere, the path of a single ray, or a bundle of rays, was analyzed using geometrical optics. These 120 Monte Carlo simulations have proved capable of assessing beam wander, in particular the (Range) 3 dependence of mean-squared beam wander, and in estimating phase shifts between the rays as they propagate through turbulence and in calculating the aperture averaging effect. An effective C n 2 was also determined by correlating beam wander behavior with the path length. The Monte Carlo simulations were all compared and show good agreement with the predictions of wave theory. This is the first report to date to simulate the aperture averaging factor using geometrical ray analysis. The aperture averaging factor was modeled in Section 5.3.3 for Gaussian as well as uniformly distributed rays. The results assessed the reduction in scintillation with increasing aperture diameter. In addition, for the uniformly distributed rays, the aperture averaging factor F showed the sharper initial decline followed by a longer tail for higher turbulence strengths as predicted by the theoretical as well as experimental results. A future direction for this work would be to develop an experimental setup that is portable to be taken outdoors closer to the ground?s surface for higher temperature fluctuations and hence strong turbulence measurements. The experimental setup presented was about 12 meters above the ground which led to only weak and intermediate turbulence results. There is a need for more strong turbulence empirical results to compare with the existing strong turbulence theory and help in the development of new theories that fits throughout all of the turbulence levels from weak to strong. For the geometrical Monte Carlo simulations, the simulation runs were limited to 20% spheres? coverage, 1000 Rays and 100 Runs. This limitation was due to the 121 long run time which was about 40 minutes for each choice of the bubbles? standard deviation of their refractive index fluctuations. Such a program could be further optimized to reduce run time, as well as computers with faster processing time or super computers could be used to increase the number of simulation runs. This would allow the aperture averaging curves to be smoother giving better results as well as allows the study of more turbulence levels through a larger variation in the simulation parameters. 122 REFERENCES [1] J. W. Strobehn (ed.). Topics in Applied Physics, Vol. 25, Laser Beam Propagation in the Atmosphere. Springer-Verlag: New York, 1978. [2] D. L. Fried, ?Aperture Averaging of Scintillation.? J. Opt. Soc. Am., Vol. 57, No. 2, February 1967, p. 169-75. [3] James H. Churnside, ?Aperture averaging of optical scintillations in the turbulent atmosphere.? Appl. Opt., Vol. 30, No. 15, 20 May 1991, p. 1982-94. [4] Larry C. Andrews and Ronald L. Phillips, Laser Beam Propagation through Random Media, SPIE Optical Engineering Press: Bellingham, Washington, 1998. [5] Albert D. Wheelon, Electromagnetic Scintillation, Vol. 2, Weak Scattering. Cambridge University Press: Cambridge, 2003. [6] James H. Churnside, ?Aperture-Averaging Factor for Optical Propagation Through the Turbulent Atmosphere.? NOAA Technical Memorandum ERL WPL-188, November 1990. (available from the National Technical Information Service ). [7] Larry C. Andrews, Ronald L. Phillips, and Cynthia Y. Hopen, ?Aperture averaging of optical scintillations: power fluctuations and the temporal spectrum,? Wave Random Media Vol. 10, 2000, p. 53-70. 123 [8] Larry C. Andrews, Ronald L. Phillips, and C. Y. Hopen, Laser Beam Scintillation with Applications, SPIE Optical Engineering Press, Bellingham, Washington, 2001. [9] Stuart D. Milner and Christopher C. Davis, ?Hybrid Free-Space Optical/RF Networks for Tactical Operations,? in Proc. IEEE, MILCOM 2004 Conference (IEEE Number: 04CH37621C, ISBN: 0-7803-8848-8). [10] H. Willebrand and B. S. Ghuman, Free-Space Optics: Enabling Optical Connectivity in Today?s Networks, SAMS, Indianapolis, 2002. [11] Akira Ishimaru, Wave Propagation and Scattering in Random Media, IEEE Press and Oxford University Press, 1997. [12] V. I. Tatarskii, The Effects of the Turbulent Atmosphere on Wave Propagation (trans. for NOAA by Israel Program for Scientific Translations, Jerusalem, 1971). [13] F. H. Champagne, C. A. Friehe, J. C. LaRue, and J. C. Wyngaard, ?Flux measurements, flux-estimation techniques, and fine-scale turbulence measurements in the unstable surface layer over land,? J. Atmosp. Sci., Vol. 34, 1977, p. 515-530. [14] R. M. Williams and C. A. Paulson, ?Microscale temperature and velocity spectra in the atmospheric boundary layer boundary layer,? J. Fluid Mech., Vol 83, 1977, 547-567. [15] Hill, ?Models of the scalar spectrum for turbulent advection,? J. Fluid Mech., Vol. 88, 1978, p. 541-662. 124 [16] James H. Churnside, ?A spectrum of refractive-index turbulence in the turbulent atmosphere,? J. Mod. Opt., Vol. 37, 1990, p. 13-16. [17] Rod G. Frehlich, ?Laser scintillation measurements of the temperature spectrum in the atmospheric surface layer,? J. Atmos. Sci., Vol. 49, 1992, p. 1494-1509. [18] Christopher C. Davis, ?Correlation Functions, Aperture Averaging, And Related Issues For Line-Of-Sight Optical Links Through Weak and Strong Scattering,? A Tutorial Discussion, Not Published to date. [19] R.W. Fante, ?Electromagnetic beam propagation in turbulent media,? Proc. IEEE, Vol. 63, 1975, p. 1669-1692. [20] R.S. Lawrence and J.W. Strohbehn, ?A survey of clear-air propagation effects relevant to optical communications,? Proc. IEEE, Vol. 58, 1970, p. 1523- 1545. [21] L.A. Chernov, Wave Propagation in a Random Medium, McGraw-Hill, New York, 1960. Reprinted by Dover, New York, 1967. [22] D. L. Fried ?Optical Resolution Through a Randomly Inhomogeneous Medium for Very Long and Very Short Exposures.? J. Opt. Soc. Am., Vol. 56, No. 10, October 1996, p. 1372-9. [23] R. L. Fante, ?Electric field spectrum and intensity covariance of a wave in a random medium,? Radio Science, Vol. 10, 1975, p. 77-85. [24] R. L. Fante, ?Irradiance scintillations: comparison of theory with experiment,? J. Opt. Soc. Am., Vol. 65, 1975, p. 548-550. 125 [25] I. S. Gradshteyn and I. W. Ryzhik, Tables of Integrals, Series and Products, Academic Press, New York, 1965. [26] Larry C. Andrews, Ronald L. Phillips, Cynthia Y. Hopen, and M. A. Al- Habash, ?Theory of optical scintillations,? J. Opt. Soc. Am. A, Vol. 16, 1999, p. 1417-1429. [27] V. I. Tatarski, Wave Propagation in a Turbulent Medium, McGraw-Hill, New York, 1961. Reprinted by Dover, New York, 1967. [28] Robert S. Lawrence and John W. Strobehn. ?A Survey of Clear-Air Propagation Effects Relevant to Optical Communications.? Proc. IEEE, Vol. 58, No. 10, October 1970, p. 1523-44. [29] Linda M. Wasiczko ?Techniques to Mitigate the Effects of Atmospheric Turbulence on Free Space Optical Communication Links,? Ph.D. Thesis, University of Maryland, 2004. [30] C. I. Moore, H. R. Burris, M. R. Suite, M. F. Stell, M. J. Vilcheck, M. A. Davis, R. Mahon, W. S. Rabinovich, G. C. Gilbreath, E. Oh, W. J. Scharpf, and A. E. Reed, ?Spatial intensity correlation and aperture averaging measurements in a 20 mile retro-reflected lasercom link,? Proc. SPIE, Vol. 5160, Free-Space Laser Communication and Active Laser Illumination III, 2004, p. 474-482. [31] H. T. Yura and W. G. McKinley, ?Optical scintillation statistics for IR ground-to-space laser communication systems,? Appl. Opt, Vol. 22, 1983, p. 3353-3358. 126 [32] Christopher C. Davis and Igor I. Smolyaninov, ?The effect of atmospheric turbulence on bit-error-rate in an on-off-keyed optical wireless system,? Proc. SPIE, Vol. 4489, Free-Space Laser Communication and Laser Imaging, 2002, p. 126-137. [33] Heba Yuksel and Christopher C. Davis. ?A Geometrical Optics Approach for Modeling Atmospheric Turbulence.? Proc. SPIE, Vol. 5891, Atmospheric Optical Modeling, Measurement, and Simulation Conference, August 2005. [34] Walid A. Atia, ?Theoretical and Experimental Investigations of Laser Beam Propagation in Turbulent Media,? M.S. Thesis, University of Maryland, 1993. [35] S. F. Clifford, ?The Classical Theory of Wave Propagation in a Turbulent Medium.? In: Laser Beam Propagation in the Atmosphere. J. W. Strohbehn (Editor), Springer-Verlag, Berlin, 1978, p. 9-41. [36] R. W. Lee, and J. C. Harp. ?Weak Scattering in Random Media, with Applications to Remote Probing.? Proc. IEEE, Vol. 57, 1960, p. 375-406. [37] J. W. Strohbehn, ?Modern Theories in the Propagation of Optical Waves in a Turbulent Medium.? In: Laser Beam Propagation in the Atmosphere. J. W. Strohbehn (Editor), Springer-Verlag, Berlin, 1978, p. 45-104. [38] B. T. Morgan, Elements of Simulation, Cambridge: The University Press, 1984. 127