ABSTRACT Title of Dissertation: LASER BEAM CONTROL, COMBINING, AND PROPAGATION IN ATMOSPHERIC TURBULENCE William R. Nelson, Doctor of Philosophy, 2016 Dissertation directed by: Professor Christopher C. Davis Department of Electrical and Computer Engineering This dissertation is concerned with the control, combining, and propagation of laser beams through a turbulent atmosphere. In the first part we consider adaptive optics: the process of controlling the beam based on information of the current state of the turbulence. If the target is cooperative and provides a coherent return beam, the phase measured near the beam transmitter and adaptive optics can, in principle, correct these fluctuations. However, for many applications, the target is uncooperative. In this case, we show that an incoherent return from the target can be used instead. Using the principle of reciprocity, we derive a novel relation between the field at the target and the scattered field at a detector. We then demonstrate through simulation that an adaptive optics system can utilize this relation to focus a beam through atmospheric turbulence onto a rough surface. In the second part we consider beam combining. To achieve the power levels needed for directed energy applications it is necessary to combine a large number of lasers into a single beam. The large linewidths inherent in high-power fiber and slab lasers cause random phase and intensity fluctuations occurring on sub-nanosecond time scales. We demonstrate that this presents a challenging problem when attempting to phase-lock high-power lasers. Furthermore, we show that even if instruments are developed that can precisely control the phase of high-power lasers; coherent combining is problematic for DE applications. The dephasing effects of atmospheric turbulence typically encountered in DE applications will degrade the coherent properties of the beam before it reaches the target. Finally, we investigate the propagation of Bessel and Airy beams through atmospheric turbulence. It has been proposed that these quasi-non-diffracting beams could be resistant to the effects of atmospheric turbulence. However, we find that atmospheric turbulence disrupts the quasi-non-diffracting nature of Bessel and Airy beams when the transverse coherence length nears the initial aperture diameter or diagonal respectively. The turbulence induced transverse phase distortion limits the effectiveness of Bessel and Airy beams for applications requiring propagation over long distances in the turbulent atmosphere. LASER BEAM CONTROL AND COMBINING IN ATMOSPHERIC TURBULENCE by William R. Nelson 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 2016 Advisory Committee: Professor Christopher C. Davis, Chair Professor Phillip Sprangle Professor Julius Goldhar Professor Jeremy Munday Professor Peter Bernard © Copyright by William R. Nelson 2016 ii Publications W. Nelson, J. P. Palastro, C. Wu, and C. C. Davis, “Using an incoherent reflection to adaptively focus through atmospheric turbulence”, Optics Letters (2016) W. Nelson, P. Sprangle, and C. C. Davis, “Atmospheric Propagation and Combining of High-Power Lasers”, Applied Optics, 55 (11), (2016) W. Nelson, J. P. Palastro, C. Wu, and C. C. Davis, “Enhanced backscatter of optical beams reflected in turbulent air”, J. Opt. Soc. Am. A, 32(7), 1371-1378 (2015) W. Nelson, C. Wu, C. C. Davis, "Determining beam properties at an inaccessible plane using the reciprocity of atmospheric turbulence", Proceedings of SPIE Vol. 9614, 96140E (2015) W. Nelson, J. P. Palastro, C. C. Davis and P. Sprangle, “Propagation of Bessel and Airy beams through atmospheric turbulence”, J. Opt. Soc. Am. A, 31(3), 603-609 (2014) W. Nelson, J. P. Palastro, C. Wu, and C. C. Davis., "Enhanced backscatter of optical beams reflected in atmospheric turbulence", Proceedings of SPIE Vol. 9224, 922411 (2014) C. Wu, W. Nelson, and C. C. Davis, "3D geometric modeling and simulation of laser propagation through turbulence with plenoptic functions", Proceedings of SPIE Vol. 9224, 92240O (2014) C. Wu, J. Ko, W. Nelson, C. C. Davis, "Phase and amplitude wave front sensing and reconstruction with a modified plenoptic camera", Proceedings of SPIE Vol. 9224, 92240G (2014) C. Wu, W. Nelson, J. Ko, C. C. Davis, "Experimental results on the enhanced backscatter phenomenon and its dynamics", Proceedings of SPIE Vol. 9224, 922412 (2014) iii Acknowledgements The work contained within this dissertation is the result of collaboration between many people. I thank my advisor, Professor Chris Davis, for his insight, guidance, and encouragement throughout my studies. He is a model on how to advise students. Additionally, I was fortunate enough to have a second advisor, Professor Phil Sprangle. His immense knowledge of physics along with passion and excitement for physics made working with him an instrumental experience. I also owe great thanks to Dr. John Palastro for his immeasurable help throughout graduate school. Not only did I learn many scientific concepts from John, but he also helped me to greatly improve my scientific writing skills. I would also like to thank Professor Larry Andrews and Professor Ron Phillips at the University of Central Florida for fruitful discussions. I am also grateful for the help and advice from the Maryland Optics Group: Chensheng Wu, Johnathan Ko, Dr. John Rzasa, and Thomas Shen. Additionally, I would like to thank the committee members: Professor Julius Goldhar, Professor Jeremy Munday, and Professor Peter Bernard for taking the time to review my dissertation. This work was supported by funding from JTO through ONR under contract number N000141211029. iv Table of Contents Publications ................................................................................................................... ii Acknowledgements ...................................................................................................... iii List of Figures .............................................................................................................. vi Definitions and abbreviations ...................................................................................... xi Chapter 1 : Overview .................................................................................................... 1 1.1: Introduction ........................................................................................................ 1 1.2: Historical Background ....................................................................................... 3 1.3: Mathematical Background ................................................................................. 6 Wave propagation ................................................................................................. 6 Atmospheric turbulence theory ............................................................................. 9 Phase screen approximation ................................................................................ 12 Linear systems and the superposition integral .................................................... 15 1.4: Motivation, overview, and original contributions ........................................... 17 Adaptive optics ................................................................................................... 17 Beam Combining ................................................................................................ 19 Quasi-non-diffracting beams .............................................................................. 20 Chapter 2 : Reciprocity and Adaptive Optics ............................................................. 23 2.1: Overview .......................................................................................................... 23 2.2: Derivation of the reciprocity relation............................................................... 25 2.3: Random walk in the complex plane ................................................................. 28 2.4: Numerical verification of reciprocity relation ................................................. 30 2.5: Adaptively focusing through turbulence.......................................................... 32 2.6: Discussion on reciprocity and adaptive optics ................................................. 36 Chapter 3 : Enhanced Backscatter .............................................................................. 38 3.1: Background ...................................................................................................... 38 3.2: Lab-scale turbulence ........................................................................................ 42 Experimental arrangement .................................................................................. 42 Simulation description ........................................................................................ 44 Retroreflector ...................................................................................................... 45 Diffuse reflector .................................................................................................. 47 3.3: Directed energy scale ....................................................................................... 50 Transition to DE length scales ............................................................................ 50 Directed energy parameters ................................................................................ 53 Retroreflector ...................................................................................................... 54 Diffuse reflector .................................................................................................. 56 3.4: Tilt-Shift method.............................................................................................. 58 3.5: Discussion on enhanced backscatter ................................................................ 61 Chapter 4 : Beam Combining .................................................................................... 62 v 4.1: Overview of beam combining .......................................................................... 62 4.2: Beam Director Model ...................................................................................... 64 Geometry of the beam director ........................................................................... 64 Large Linewidth Source ...................................................................................... 65 Random intensity and phase fluctuations ........................................................... 68 Coherently combining tiles ................................................................................. 72 4.3: Propagation in turbulent atmosphere ............................................................... 75 Vacuum ............................................................................................................... 76 Moderate turbulence ........................................................................................... 79 Strong turbulence ................................................................................................ 82 4.4: Discussion on beam combining ....................................................................... 85 Chapter 5 : Bessel and Airy Beams ............................................................................ 88 5.1: Introduction ...................................................................................................... 88 5.2: Gaussian beam ................................................................................................. 89 5.3: Bessel beam ..................................................................................................... 92 Fixed aperture comparison .................................................................................. 97 5.4: Airy beam ........................................................................................................ 99 5.5: Discussion on quasi-non-diffracting beams ................................................... 103 Chapter 6 : Summary and future work ..................................................................... 105 6.1: Summary ........................................................................................................ 105 6.2: Future work .................................................................................................... 107 References ................................................................................................................. 108 vi List of Figures Figure 1.1: Illustration of a generic phase screen simulation. Advancing the beam through each segment of length z involves three steps. First the beam is propagated through vacuum a distance / 2z . Then a random phase screen is applied to the field. Finally the beam is propagated through vacuum the remaining distance / 2z . The process is repeated for each segment. ......................................................................... 14 Figure 2.1: (a) Schematic of the propagation geometry. The transmitter and receiver plane coincide, creating a monostatic channel. The medium between the transmitter/receiver plane is arbitrary, but in this context we consider a turbulent atmosphere. A thin lens is located between the receiver and detector planes, equidistant from each plane. (b) Schematic of the phase modulator located in the transmitter/receiver plane. The phase modulator consists of square elements of length a which apply piston type phase perturbations to the transmitted and received beams. The total dimension of the phase modulator is defined as b. ...................................... 24 Figure 2.2: (a) Isotropic random walk in the complex plane with 64 steps of length 1. (b) Isotropic random walk in the complex plane with 16 steps of length 4. ............... 30 Figure 2.3: . Scatter plot of the magnitude (a) and phase (b) depicting the numerical equivalent of the relation described by Eq. (2.7). The vertical and horizontal axes correspond to the left and right hand sides of Eq. (2.7) respectively. Since the values are complex, there are two scatterplots to demonstrate that the simulation in good agreement with the theory. .......................................................................................... 31 Figure 2.4: Top row: Initial (a) and final (b) intensity in the instance when SPGD was most effective. The metric increased from 304J  to 1,280J  bewtween (a) and (b). The intensity profile (a) before the SPGD algorithm is distorted by atmospheric turbulence. After 100 iterations of SPGD the intensity profile (b) displays a localized region of high intensity. .............................................................................................. 35 Figure 2.5: (a) Histogram of maximum intensity over 300 independent channels comparing the intensity profiles at iteration 1 (red) and iteration 100 (blue) of the SPGD algorithm. The average maximum intensity during iteration 1 is 1.9 with a standard deviation of .36. The average maximum intensity during iteration 100 is 3.5 with a standard deviation of 1.2. (b) Histogram of PIB contained within a 1 cm radius bucket centered on the location of maximum intensity. The average PIB during iteration 1 is .027 with a standard deviation of .005. The average PIB during iteration 100 is .045 with a standard deviation of.015. ............................................................. 36 Figure 3.1: Diagram of a pair of reciprocal rays taking the same path in opposite directions. .................................................................................................................... 41 vii Figure 3.2: Schematic of the experimental arrangement for the (a) monostatic and (b) bistatic channels. Turbulence is generated by two 61 cm hotplates for each channel. 43 Figure 3.3: Average monostatic (first column), bistatic (second column) and x-slice (third column) intensity profiles from a retroreflector target obtained from the experiment (first row) and simulation (second row). The experimental and numerical results are in good agreement. The localized region of high average intensity is observed in the monostatic channel and is increased by approximately a factor of two in comparison to the bistatic channel. ......................................................................... 46 Figure 3.4: Average monostatic (first column), bistatic (second column) and x-slice (third column) intensity profiles from a diffuse reflector target obtained from the experiment (first row) and simulation (second row). In the experimental arrangement the target is a flat, white poster board. The EBS peak is approximately twice as wide in the experiment compared to the simulation results. This is likely due to the material properties associated with the target [61]. .................................................................. 48 Figure 3.5: Comparison of lab and DE simulations for a retroreflector (a) and diffuse reflector (b). With the appropriate scaling, we are able to find equivalent large-scale parameters which correspond to the lab-scale experiment. ........................................ 52 Figure 3.6: Enhancement factor as a function of 2 nC (a) and number of independent runs averaged (b) for a retroreflector target. The trend of enhancement factor versus 2 nC is in good agreement with early theoretical predictions. Plot (b) demonstrates that many runs must be averaged over to observe the enhancement factor of 2. .............. 54 Figure 3.7: Enhancement factor as a function of 2 nC (a) and number of independent runs averaged (b) for a diffuse reflector target. The trend of enhancement factor versus 2 nC is in good agreement with early theoretical predictions. Plot (b) demonstrates that many runs must be averaged over to observe the enhancement factor of 2. ................................................................................................................... 57 Figure 3.8: Average monostatic (a) bistatic (b) and x-slice (c) intensity profiles using the tilt shift method in static turbulence. The region of localized high intensity is observed from a single instance of turbulence. This allows observation of EBS without performing a time average. ............................................................................ 60 Figure 4.1: Beam director geometry. We define the side length of each tile as a and the distance between the centers of adjacent tiles as b where b a . Each tile is indexed by the indices l and m defined on the set  1,0,1 . ...................................... 65 Figure 4.2: Comparison of a Gaussian and Lorentzian power spectrum when the linewidths are matched, i.e. G L       . The Gaussian power spectrum viii contains 98% of the total power within  of 0 . The Lorentzian power spectrum contains 70% of the total power within  of 0 . ................................................... 69 Figure 4.3: Intensity (a) and phase (b) fluctuations for a single tile with a Gaussian power spectrum. The intensity has random fluctuations when observed on time scales comparable to the coherence time............................................................................... 70 Figure 4.4: Intensity (a) and phase (b) fluctuations for a single tile with a Lorentzian power spectrum. The intensity has random fluctuations when observed on time scales comparable to the coherence time............................................................................... 71 Figure 4.5: The RMS phase difference between two tiles with a Gaussian power spectrum (blue) and Lorentzian power spectrum (red) as a function of instrument bandwidth .................................................................................................................... 74 Figure 4.6: Intensity profiles for propagation through vacuum at 5 kmz  for the incoherently combined beam (a), coherently combined beam (b), and phase matched monochromatic beam (c). Each tile of the incoherently combined beam diffracts independently resulting in significant spreading. The coherently combined beam has an initial RMS phase difference 2 6RMS   between any two tiles. The far-field intensity profile of the coherently combined beam resembles that of the monochromatic beam. ................................................................................................. 77 Figure 4.7: Vacuum intensity profiles along x-axis at for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). Intensity is normalized by the on-axis intensity of the monochromatic beam. The on axis intensity of the incoherently and coherently combined beams are approximately 11% and 71% of the on-axis intensity of the monochromatic beam respectively. ................................................................................................................ 78 Figure 4.8: Power in the bucket (PIB) at 5 kmz  through vacuum as a function of bucket radius for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). The coherently combined beam and monochromatic beam contain significantly more power near the axis due to the distinct central lobe. The contained power within a radius of 20 cm is approximately equal for all beams. ............................................................................. 79 Figure 4.9: Intensity profiles for propagation through moderate turbulence, 2 14 2 310 mnC   , at 5 kmz  for the incoherently combined beam (a), a coherently combined beam (b), and phase matched monochromatic beam (c). The coherently combined beam and monochromatic beam have smaller spot size and higher on-axis intensity than the incoherently combined beam.......................................................... 80 Figure 4.10: Average intensity profiles in moderate turbulence, 2 14 2 310 mnC   , along x-axis at 5 kmz  for the incoherently combined beam (blue), coherently ix combined beam (green), and phase matched monochromatic beam (red). Intensity is normalized by the on-axis intensity of a phase-matched, monochromatic beam propagated through vacuum. For reference, the intensity profile of a phase-matched, monochromatic beam propagated through vacuum is denoted by the black curve. The on-axis intensity of the coherently combined beam is approximately twice the on-axis intensity of the incoherently combined beam. ............................................................ 81 Figure 4.11: Power in the bucket (PIB) at 5 kmz  through moderate turbulence, 2 14 2 310 mnC   , as a function of bucket radius for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). For reference, the PIB of the monochromatic, phase matched beam propagated through vacuum is denoted by the dashed black curve. The coherently combined beam and monochromatic beam contain more power near the axis than the incoherently combined beam, however the difference is not as significant as in vacuum. ....................................................................................................................... 82 Figure 4.12: Intensity profiles for propagation through strong turbulence, 2 13 2 310 mnC   , at 5 kmz  for the incoherently combined beam (a), coherently combined beam (b), and phase matched monochromatic beam (c). The intensity profiles are nearly indistinguishable. .......................................................................... 83 Figure 4.13: Intensity profiles in strong turbulence, 2 13 2 310 mnC   , along x-axis at 5 kmz  for the incoherently combined beam (blue), coherently combined beam (green), monochromatic beam (red). Intensity is normalized by the on-axis intensity of a monochromatic beam in vacuum. The intensity profiles are nearly indistinguishable. ........................................................................................................ 84 Figure 4.14: Power in the bucket (PIB) at 5 kmz  through strong turbulence, 2 13 2 310 mnC   , as a function of bucket radius for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). All beams deliver approximately the same power. ........................................... 85 Figure 5.1: Transverse profiles of a Gaussian beam with initial spot size 0 4.5w cm after 6.4 km of propagation in (a) vacuum, (b) turbulence with 2 15 2/31 10nC m    , and (c) averaged over 100 instances through turbulence. ........................................... 90 Figure 5.2: Simulation results for a Gaussian beam. (a) and (b) Comparison of the on- axis intensity in turbulence, red, and vacuum, blue, as a function of propagation distance for initial spot sizes of 1.8 cm and 4.5 cm respectively. (c) Ratio of the RMS radius at one Rayleigh length to the initial spot size as a function of initial spot size. (d) Normalized on axis intensity at one Rayleigh length as a function of initial spot size. In (a-d) the ratio of the Fried parameter to beam diameter, green, is plotted for reference. ..................................................................................................................... 92 x Figure 5.3: Transverse profiles of a 15 ring Bessel beam hard-apertured at a radius of 23 cm after 6.4 km of propagation in (a) vacuum, (b) turbulence with 2 15 2/31 10nC m    , and (c) averaged over 100 instances through turbulence. .......... 94 Figure 5.4: Simulation results for a Bessel beam. (a) and (b) Comparison of the on- axis intensity in turbulence, red, and vacuum, blue, as a function of propagation distance for initial aperture radii of 6.65 cm and 23 cm respectively. (c) Ratio of the RMS radius at one diffraction length, BL , to the initial aperture radius as a function of initial aperture radius. (d) Normalized on axis intensity at one diffraction length as a function of initial aperture radius. In (a-d) the ratio of the Fried parameter to beam diameter, green, is plotted for reference. .................................................................... 96 Figure 5.5: The power delivered to a circular aperture of radius 15 cm (top) and on axis intensity (bottom) as a function of rings in a Bessel beam initially apertured at a radius of 15 cm after propagation distances of 1.6 km, blue, 4 km, red, and 6.4 km, green. The dashed and solid lines are results from propagation in vacuum and turbulence with 2 15 2/31 10nC m    respectively. ...................................................... 98 Figure 5.6: Transverse profiles of a 15 zero Airy beam with 2 22 31.3d cm cm   after 6.4 km of propagation in (a) vacuum, (b) turbulence with 2 15 2/31 10nC m    , and (c) averaged over 100 instances through turbulence. 101 Figure 5.7: Simulation results for an Airy beam. (a) and (b) Comparison of the on- axis intensity in turbulence, red, and vacuum, blue, as a function of propagation distance for an aperture diagonal of 11.3 cm and 31.1 cm respectively. (c) Ratio of the RMS radius at one diffraction length, AL , to the aperture diagonal as a function of aperture diagonal. (d) Normalized on axis intensity at one diffraction as a function of aperture diagonal. In (a-d) the ratio of the Fried parameter to beam diagonal, green, is plotted for reference. ................................................................................................. 102 xi Definitions and abbreviations bistatic channel A channel in which the beam propagates through different media on the paths from the transmitter to the target and the target to the receiver. EBS Enhanced backscatter. A phenomenon characterized by a region of localized high intensity after a double pass through a random medium enhancement factor Ratio of average intensities from a monostatic channel and bistatic channel at the enhanced backscatter location FWHM full width half maximum monostatic channel A channel in which the beam propagates through the same medium on the paths from the transmitter to the target and the target to the receiver PIB power in the bucket RMS Root mean squared SPGD Stochastic parallel gradient descent. An algorithm designed to optimize a metric by applying random perturbations to all controls at once and calculating the change in the metric TSM Tilt shift method. A novel algorithm to expedite the process of observing enhanced backscatter. It replaces the standard temporal average with a spatial average. 1 Chapter 1 : Overview 1.1: Introduction Propagation of waves through random media has long been an active area of research. It was first investigated in the context of propagation of starlight, sound, and microwaves through the atmosphere and oceans. With the invention of the laser in 1960, the number of applications involving propagation of waves through random media significantly increased. Some of these new applications include directed energy, remote sensing, free space optical communications, and power beaming. The common problem faced by all these application and the topic of this dissertation, is propagating a laser beam long distances through the atmosphere while maintaining beam quality. Phase distortions acquired from turbulent fluctuations in the refractive index of the atmosphere modify the propagation of laser beams [1, 2, 3, 4, 5, 6, 7]. Although these fluctuations in the refractive index occur on a scale of 10 -6 , the distortions become profound when considering multiple kilometer propagation ranges. In general, the modifications to the beam can be described as three effects: wander, spreading, and scintillation, each of which is detrimental to applications requiring long range propagation of laser beams [4, 7]. For instance, beam wander, or deflections in the beam’s centroid, can result in the beam missing a target entirely. Beam spreading, expansion of the beam beyond vacuum diffraction, distributes the beam power over a larger area, reducing the intensity, delivered power, and efficiency. And scintillation, 2 fluctuations in the beam’s intensity, can result in image distortion or higher bit-error- rates in optical communication systems [8]. A laser beam can be controlled in a manner which reduces the distortions caused by atmospheric turbulence. In theory, this control could take place with or without information about the current state of the atmosphere. Actively controlling the radiation based on the current state of the propagation medium is termed adaptive optics. One common implementation of adaptive optics uses a beacon beam to gain information about channel (ideally this information would be in the form of the point spread function of the channel) and a deformable mirror to apply corrections to the light in a real-time, closed-loop fashion. In contrast, a proposed method which does not use information on the current state of the turbulence is the propagation of quasi- non-diffracting beams. This method consists of generating Bessel or Airy beams and relying on the ‘self-healing’ properties of the beam to resist the negative effects of atmospheric turbulence. As we will demonstrate, this is not an effective method for propagating a beam through turbulence. Some applications of laser beam propagation through the atmosphere, such as directed energy, require high powers, which cannot be achieved by a single fiber or slab laser. In this case, multiple laser beams must be combined. These beams can be combined coherently or incoherently. Coherent beam combining involves precisely controlling the phases of each beam so that they fluctuate synchronously. The motivation behind this phase synchronization is that it allows the combination of many beams to propagate as a single large beam which then reduces the diffractive spreading angle and in theory, the beam can be focused to a smaller spot. However, as 3 we will demonstrate, the benefits of coherent combining diminish in moderate to strong turbulence. A much simpler approach with comparable results in moderate and strong turbulence is incoherent beam combining. This introductory chapter is intended to give some historical and mathematical background of laser beam propagation through the atmosphere. Additionally, there are descriptions of the specific topics covered, motivations for the research, and brief summaries of the findings and results. 1.2: Historical Background Before the 1950’s astronomers widely believed that the limitations imposed on astronomical observations by atmospheric turbulence could not be overcome. Then in 1953, Babcock proposed the first adaptive optics system [9]. It utilized an Eidophor, a type of television projector, as the deformable optical element. Eidophors consist of a thin oil layer on a mirror to which a modulated electric charge is transferred. It was effectively the first form of a deformable mirror. In 1957, Linnick described how a beacon or guide star could be used to probe the atmosphere and control a deformable optical element [10]. These two ideas prompted further study of waves propagating through random media. Although most early studies were conducted for astronomical imaging, the findings are still relevant to horizontal propagation of laser beams. Of special importance was the work carried out by Fried, Greenwood, and Primmerman among others. Fried introduced the notion of the transverse coherence length, or Fried parameter, along with concepts of short and long exposure times [11]. The Fried parameter is the diameter of a circle over which the expected RMS phase difference is 1. In the following section, the Fried parameter and its implications will be described 4 in greater detail. The studies of Greenwood in 1977 determined the minimum requirement on the speed response of an adaptive optics system feedback loop based on atmospheric conditions [12]. In general, an adaptive optics system should operate on a scale of at least kHz, however there are many factors which influence this. Some of these factors will be discussed in later sections. Prior to Primmerman’s work in 1991, astronomers would use a bright star near the object they wished to image as a guide star for the deformable mirror. Clearly this was not ideal, as astronomers were restricted to viewing objects that had a sufficiently bright star within the isoplanatic angle, or the angular distance over which the turbulence is nearly unchanged [1]. Primmerman presented the first experimental results on using a laser beam as an artificial guide star [13]. The idea is based on propagating a laser high into the atmosphere and using the backscattered light as a guide star to control the adaptive optics system. The incorporation of lasers into adaptive optics system is not limited to providing an artificial guide star for astronomical observations. Adaptive optics systems started being researched for inclusion in free space optical communication systems and directed energy systems. These applications differ from astronomical imaging in that the radiation is coherent and the propagation path is mainly in the horizontal direction. During the 1980’s and 1990’s there was a strong desire to have a laser defense system which could disable incoming ballistic missiles. A laser defense system faces two difficulties: achieving power levels high enough to damage a target and the degradation of the beam due to atmospheric turbulence. These two are related, as being able to compensate for atmospheric turbulence can reduce the 5 minimum power requirements of the laser. The US Air force initiated the Airborne Laser Project in 1996 which led to the development of a laser with enough power to accomplish such a task. To achieve high powers, the developers used a chemical oxygen iodine laser which was too bulky and required too much fuel to be practical. It was so large it could only be carried on a Boeing 747 and there was little space left for laser fuel [14]. The solution to the size problem was solved by recent advances in fiber lasers by companies such as IPG photonics and Nufern [7]. However, to achieve the power levels needed it is necessary to combine a large number of lasers into a single beam. The beams can be combined coherently or incoherently. In chapter 5 we will discuss and compare these beam combining strategies in detail. Once the lasers are combined into a single beam, adaptive optics can be used to increase the efficiency and extend the range of the beam. However, contrary to astronomical imaging, it is not clear what to use as a guide star because using the backscatter from the upper atmosphere is not appropriate for horizontal propagation of the beam. In bidirectional free-space optics communication systems, a beacon laser can be propagated between nodes to gain information about the turbulence and control the deformable mirror. This does not apply to instances where the target is uncooperative and it remains unclear how to control the deformable mirror. We will address this situation in detail in chapter 2 where we simulate a novel design for an adaptive optics system in which the target is uncooperative. 6 1.3: Mathematical Background Wave propagation We begin with Maxwell’s equations in MKS units with vacuum magnetic permeability and absent of free charges: 0 t       H E (1.1) t      E H (1.2) 0 E (1.3) 0 0 H , (1.4) where E is the vector electric field, H is the vector magnetic field, 0 is the vacuum magnetic permeability, and  is the electric permittivity of the medium. We derive the wave equation by taking the cross product of both sides of Eq. (1.1) and applying the vector identity 2( ) ( )    E E E to the left side: 2 0 ( ) ( ) t          H E E . (1.5) Substituting Eq. (1.2) into the right side of Eq. (1.5) provides 2 2 0 2 ( ) t          E E E . (1.6) We consider the medium to be inhomogeneous so that the electric permittivity is a function of space. Using the vector identity ( )      E E E we can rewrite Eq. (1.3) as       E E . (1.7) Substituting Eq. (1.7) into Eq. (1.6) and moving all terms to the left side gives the result 7 2 2 0 2 0 t               E E E . (1.8) At this point, it is convenient to define the velocity of light in vacuum, which is 0 0 1 c    , (1.9) and the refractive index 1 2 0 n          , (1.10) where 0 is the vacuum electric permittivity. We can express the refractive index of air as 1n n  where 1n . This allows us to make the approximation 2 1 2n n  and consequently 0(1 2 )n    . We can further express n as the sum of two components, ( )Tn n n    r , where n accounts for the homogenous index based on the chemical compounds of air and ( )Tn r accounts for the inhomogeneous fluctuations due to turbulence with ( ) 0Tn r . The average, denoted by , is performed over an ensemble of statistically independent instances of the index fluctuations. Later, we will discuss the term ( )Tn r in greater detail. We can now approximate Eq. (1.8) as   2 2 2 2 1 2 ( ) 1 2 2 ( ) 0T Tn n n c t                 E E E r r . (1.11) We express E as a plane wave propagating in the z direction modulating a slowly varying envelope E :  1 0 02( , , ) ( , , )exp . .x y z x y z i k z t c c    E E (1.12) where 0 02 /k   and 0 0ck  are the carrier wave number and central frequency of the radiation respectively. We divide the Laplacian into transverse and forward components such that 2 2 2 2z    and use Eq. (1.12) to rewrite Eq. (1.11) as 8     2 2 2 0 02 2 2 0 02 2 0 0 2 1 1 2 2 ( ) 2 2 ( ) exp 0 T T ik k z z n n i c t t n i k z t                                              E r E E r (1.13) We can ignore the time derivatives if we assume that the envelope is in steady state. We can also ignore the second derivative in z if we assume that the envelope is slowly varying in z such that 0z k E E . This is called the paraxial approximation. Equation (1.13) is now simplified as     2 0 2 0 0 0 2 2 ( ) 2 ( ) exp 0 T T ik z k n n n i k z t                             E r E E r (1.14) The third term in Eq. (1.14) which couples the vector components can also be neglected. This is apparent when we compare the magnitudes of the third term  ( )Tn E r with the magnitude of the second term 2 0 ( )Tk n r E . From Eq. (1.12) it is clear that E E so we simply need to compare ( )Tn r with 2 0 ( )Tk n r . If we assume that the wavelength of radiation is much smaller than the characteristic length scale of the turbulent fluctuations then we have the relation 0T Tn k n  . Thus we can conclude that 2 0T Tn k n  and   2 0( ) ( )T Tn k n  E r r E . Neglecting the third term of Eq. (1.14), the relation is reduced to simply 2 2 0 02 2 ( )Tik k n n z               E r E (1.15) 9 We can further simplify Eq. (1.15) by expanding E as  ˆ exp i kzE E where k is the small change in wavenumber due to the homogenous part of the refractive index, 0k k n  . Straight forward substitution and differentiation leads to 2 2 0 0 ˆ ˆ2 2 ( )Tik k n z          E r E (1.16) where  1 02 ˆ( , , ) ( , , )exp . .x y z x y z i kz t c c    E E (1.17) and 0k k k  . The vector components of Eˆ are now clearly separable in Eq. (1.16) and we are able to write the differential equation in scalar form. Changing EE where E is the scalar electric field in a transverse direction and ˆ EE where E is the electric field envelope, we can write the scalar versions of Eq. (1.16) and Eq. (1.17) as 2 2 0 02 ( ) 2 ( ) ( )Tik E k n E z          r r r , (1.18) and 1 02 ( , ) ( )exp[ ( )] . .E t E ,t i kz t c c   r r (1.19) respectively. Atmospheric turbulence theory We now turn our attention to the turbulent refractive index fluctuations, ( )Tn r , which arise from small temperature fluctuations in the air. We note that strictly speaking the refractive index depends on the atmospheric density, but because the fluctuations are nearly isobaric, the density fluctuations are directly proportional to the temperature fluctuations. In the Kolmogorov cascade theory of turbulence, 10 temperature fluctuations are formed with large scales sizes, defined as the outer scale 0L [1]. Fluctuations occurring over a distance of approximately 0L contribute primarily to beam wander. These large fluctuations continually dissipate to smaller scales through molecular diffusion and viscosity. The cascading to smaller scales occurs because the rate of dissipation increases as the fluctuation size decreases. When the dissipation rate equals the heating rate, the cascade terminates and defines the inner scale length of the fluctuations, 0 . Fluctuations at this smaller scale contribute primarily to scintillation. The refractive index fluctuations are characterized by their covariance function, ( , ) ( ) ( )n T TB n n  r r r r and the refractive index structure function, 2( , ) [ ( ) ( )]n T TD n n   r r r r [4]. Throughout the dissertation, we consider homogeneous, isotropic, Gaussian fluctuations such that ( , ) (| |)n nB B  r r r r or ( , ) 2[ (0) (| |)]n n nD B B   r r r r fully determine the statistical properties. The Fourier transform of the covariance function, ( )n  , represents the distribution of fluctuation scale sizes. This distribution can be modeled as the Kolmogorov spectrum 2 11 3( ) 0.033n nC    , (1.20) the Von Karman spectrum 20( /2 ) 2 2 2 11/6 0 ( ) 0.033 ( ) n n e C L          , (1.21) or modified Von Karman spectrum (also called the Hill spectrum)     7 6 2 2 2 11 6 2 2 0 exp ( ) 0.033 1 1.802 0.254 l n n l l C                              (1.22) 11 where 03.3l l  , 0 02 L  , and 2 nC is the refractive index structure constant. Using Eq. (1.21), the index structure function can be shown to be 2 2/3| |n nD C r r for 0 0| | L  r r which is a form of the Kolmogorov 2/3 power law [1]. To assess the modifications to a laser beam in atmospheric turbulence, it is common to consider the solution to Eq. (1.18) in the Fraunhofer diffraction limit. This limit is applicable when considering propagation distances longer than the diffraction length or propagation to the focal point of a lens. In particular, the ensemble averaged, on-axis intensity, 21 02 (0, ) | (0, ) |I z c E z , far from the aperture can be expressed 21 2 0 2 [ ( , ) ( , )]*01(0, ) ( ,0) ( ,0) 2 2 z z c k I z E E e d d z                         r r r r r r (1.23) where 0 0 ( , ) ( , ) z Tz k n z dz     r r (1.24) The exponential in Eq. (1.23) includes the phase structure function 2(| |) [ ( , ) ( , )]sD z z         r r r r . Using the refractive index structure function, one can show that 2 2 5/3 02.91 | r r |s nD k C z   for 0 0| - | L r r [1, 4]. With Eq. (1.23) and sD , we can derive an estimate of the distance, TL , at which turbulence becomes important. The ensemble averaged intensity drops by one e-folding due to turbulence when 2sD  . Setting | r r | L    , the characteristic transverse variation in the envelope, we find 2 2 5/3 1 00.69( )T nL C k L   . Up to a numerical factor this is the same distance at which the Fried parameter 2 2 3/50 0( ) 1.67( )nr z C k z  equals L [1, 4, 12 7]. The Fried parameter, also referred to as the transverse coherence length, describes the transverse distance over which the phase fronts of the laser beam remain correlated. In other words, transverse distortions on the size of 0r and modifications to beam propagation can be expected when the Fried parameter is smaller than the beam diameter. Phase screen approximation For intervals of propagation, sz , much shorter than the diffraction length, 2 0 /dL w  (for a Gaussian beam of spotsize 0w ) . the index fluctuations result in the accumulation of a transverse phase: (r , )ˆ ˆ(r , ) (r , ) si zsE z z E z e        , where 0( , ) ( , ) sz z s T z z k n z dz        r r (1.25) We consider ( , )sz  r to be a statistically homogenous random field and can represent it in terms of the Riemann-Stieltjes integral ( , ) ( ) i sz e dv           κ r r κ (1.26) Where ( )dv κ denotes a random complex amplitude. Using the form in Eq.(1.26) we can write the correlation function of accumulated phase as ( , ) ( , ) exp[ ( )] ( ) ( )s sz z i dv dv                     r r κ r κ r κ κ (1.27) To satisfy the assumption of homogeneity, we must have relation ( ) ( ) ( ) ( )dv dv d d          κ κ κ κ κ κ κ (1.28) where ( )  κ is the power spectrum of accumulated phase. Using Eq.(1.28) we can simplify Eq.(1.27) as 13 ( ) ( , ) ( , ) ( ) i s sz z e d                 κ r r r r κ κ (1.29) Additionally, we can calculate the covariance function of accumulated phase directly from Eq.(1.25), in which case the result is 0 0 0 0 2 0( , ) ( , ) ( , ) ( , ) s sz z z z s s T T z z z z k n z n z dzdz                 r r r r . (1.30) Recalling that the covariance is the Fourier transform of the power spectral density, we can write Eq.(1.30) in terms of ( )n  : 0 0 0 0 ( ) ( )2 0 ( , ) ( , ) ( ) s s z s s z z z z i i z z z n z z z z k dz dz d e d e                           κ r r r r κ κ . (1.31) Where we have separated the transverse components and z component of the Fourier transform. Rearranging the order of integration, we now consider the result of the integral 0 0 0 0 ( ) 2 2 [1 cos( )] s s z z z z z i z z z s zz z dz dz e z          (1.32) which is obtained through direct integration. In the limit z is much greater than the outer scale of turbulence, 0sz L  , it is clear that only small z contribute and we can make the approximation 2 2 [1 cos( )] 2 ( )z s s z z z z         . (1.33) Substituting the approximation in Eq. (1.33) into Eq. (1.31) we obtain ( )2 0( , ) ( , ) 2 ( ,0) i s s s nz z k z d e                  κ r r r r κ κ . (1.34) Through comparison of Eq. (1.29) and Eq. (1.34), we can determine the power spectral density of accumulated phase fluctuations 2 0( ) 2 ( ,0)s nk z     κ κ (1.35) 14 Using Eq.(1.35) we can write the Riemann-Stieltjes integral in Eq.(1.26) as 1/2 1/2 0( , ) (2 ) [ ( ) ( )] ( ,0) i s s r i nz z k d e a ia                 r r     (1.36) where ra and ia are independent standard normal Gaussian random variables. Equation (1.36) provides the basis for the phase screen approximation in which the laser beam acquires the accumulated transverse phase at discrete axial points along the propagation path, and is propagated in vacuum between these points. Figure 1.1: Illustration of a generic phase screen simulation. Advancing the beam through each segment of length z involves three steps. First the beam is propagated through vacuum a distance / 2z . Then a random phase screen is applied to the field. Finally the beam is propagated through vacuum the remaining distance / 2z . The process is repeated for each segment. The phase screen simulation technique is used frequently throughout this dissertation. To simulate the propagation, we use a pseudo-spectral, split step algorithm illustrated in Figure 1.1. A typical axial advance, propagating the beam forward by z , involves three steps. In the first step, the transverse Fourier transform of the envelope is advanced a half step, / 2z by applying the diffraction propagator: 1 2 21 02 ( )1 2 ( , ) ( , ) x y i k k k z E z z E z e        k k , where the overbar denotes the transform 15 with respect to the transverse plane. For the second step, the phase screen is applied in coordinate space ( , )1 1 2 2 ( , ) ( , ) i z E z z E z z e            r r r . Finally the beam is advanced the second half step in the transverse Fourier domain, 1 2 21 02 ( )1 2 ( , ) ( , ) x y i k k k z E z z E z z e         k k . We note that a phase screen does not need to be applied at every advance but the interval between applications should satisfy the constraint that 0 s dL z L   . Linear systems and the superposition integral Due to the random nature of the propagation medium it is sometimes convenient to describe the process in an abstract method as a linear system. For example, we define the system  S which has the property of linearity           S ap bq aS p bS q  1 1 1 1r r r r (1.37) for all input functions p, q and all complex constants a, b. For our purpose, the inputs and outputs of the operator  S are complex-valued functions describing the electric field. For example,  2 2 1( ) ( )E S E 1r r (1.38) describes propagation of the field between the parallel planes 1r and 2r where 1 1( )E r and 2 2( )E r denote the electric field envelopes at the respective planes. Using the Dirac delta function we can decompose the field 1 1( )E r as a sum over the field at individual points: 16 1 1 1 1( ) ( ) ( )E E d r ρ r ρ ρ . (1.39) The quantity 1( )E ρ in Eq. (1.39) is now merely a weighting factor analogous to a or b in Eq.(1.37). Plugging equation Eq. (1.39) into equation Eq. (1.38) results in  2 2 1 1( ) ( ) ( )E S E d r ρ r ρ ρ . (1.40) Using the property of linearity, we can move the operator inside the integral and rewrite the relation as  2 2 1 1( ) ( ) ( )E E S d r ρ r ρ ρ . (1.41) In this form, it is convenient to define a new function     2 1;G S  r ρ r ρ providing the relation 2 2 1 2( ) ( ) ( ; )E E G d r ρ r ρ ρ . (1.42) The function  2;G r ρ is called the impulse response or the point spread function and it describes the response of the system to a point source. When boundary or initial conditions are satisfied it is also called the Green’s function. Equation (1.42) is essentially the Huygens-Fresnel principle. Propagation of the electric field envelope through a homogenous medium is described by     2 0 2; exp 2 2 ik ik G z z         2r ρ r ρ (1.43) where z is the distance between the planes 1r and 2r . For propagation through an inhomogeneous medium, the Green’s function is undefined and it is used in an 17 abstract sense. The mathematics described here provides the basis for our derivation in Chapter 2. 1.4: Motivation, overview, and original contributions In this dissertation, we investigate three topics on the propagation of optical beams through atmospheric turbulence: adaptive optics, beam combining, and quasi-non- diffracting beams. Each of these topics could be implemented in conjunction with one another; however it is best to consider each individually. Adaptive optics Adaptive optics is a method of improving an optical signal by using information about the medium through which the radiation passes. It has been shown to be successful in extending the range of beam propagation in certain instances. One adaptive optics configuration utilizes a cooperative target that provides a spatially coherent emitter return from a point on or near the target. This is essentially a beacon or guide star [15, 16]. A wavefront sensor, placed near the beam transmitter, measures the field received from this emitter. If the conjugated field is transmitted, it reproduces the emitter’s irradiance profile at the target. In the special case the target emitter is point-like, the measured field provides the turbulent channel’s point spread function, and the target can be illuminated with a desired irradiance profile in the vicinity of the emitter. Often there is no coherent emitter at or near the target, the state of the turbulent channel remains unknown, and one must rely on a reflection from the target. For normal-incidence mirror reflections, the light received at the wavefront sensor has 18 passed through the turbulent channel twice, and the phase aberration acquired in each pass cannot be explicitly determined. Scattering from a rough surface adds a further complication. These surfaces scramble the beam phase such that the reflected beam loses its spatial coherence. We consider this problem in chapter 2, where we present a design for an adaptive optics system that can use the incoherent return to correct phase aberrations and focus the beam onto the target. In particular, we use the principle of reciprocity to derive a novel relation between the field incident upon a target and the field detected by a point receiver. The derivation considers a double pass through the turbulent channel: propagation from the source to the target and back to the detector, making it especially suitable for uncooperative targets. Simulations of an adaptive optics system using the relation demonstrate enhancements in both the on-target intensity and power in the bucket. Due to the double pass and reciprocal geometry considered in chapter 2, the derived relation is related to a phenomenon called enhanced backscatter (EBS). Enhanced backscatter is observed as a retroreflected, localized region of high intensity. In chapter 3, we perform an in-depth experimental and numerical study of enhanced backscatter. Additionally, we present a novel algorithm called the ‘tilt shift method,’ which greatly reduces the time required to observe the region of high intensity. Enhanced backscatter is useful in accurately determining the optical axis of a system, and would be a critical part of implementing the adaptive optics system described in chapter 2. 19 Beam Combining Advances in solid state lasers, especially slab and fiber lasers, have made them candidates for directed-energy applications [14, 17]. To achieve the power levels needed for these applications it is necessary to combine a large number of lasers into a single beam. The combined beam is then propagated many kilometers through a turbulent atmosphere. Recently, much effort has been expended in coherent beam combining [17, 18, 19, 20]. Coherent beam combining refers to matching the phases of multiple lasers either in the transmitter plane [17, 18] or the target plane [19]. However, there are a number of important issues to be considered before a coherent combining architecture can be used for DE applications. These issues include the consideration of sources with a finite spectral linewidth and the dephasing effect of atmospheric turbulence. In chapter 4 we investigate these issues faced by coherent beam combining and discuss the physical processes associated with the combining and propagation of high-power laser beams. We compare the energy delivered to a target for the case of coherently combined and incoherently combined laser beams. Through simulation, we demonstrate that in vacuum and weak turbulence, the effectiveness of coherent combining is limited by the phase fluctuations occurring due to the broad linewidth of high-power lasers (Figure 4.7 and Figure 4.8). Furthermore, the advantages of coherent combining begin to diminish in moderate turbulence, as the turbulence induced phase distortions become dominant (Figure 4.10 and Figure 4.11). In strong turbulence and multi-km propagation ranges, we find negligible differences in the energy on the target between coherently and incoherently combined laser beams 20 (Figure 4.12 and Figure 4.13). The incoherently combined architecture is far simpler to implement and both beam combining architectures have the capability of employing adaptive optics to extend the range. Additionally, some incoherent combining methods, such as spectral beam combining, can spatially overlap all the beams so that they travel through the same turbulent fluctuations and only one adaptive optics element is needed to apply the corresponding corrections. Quasi-non-diffracting beams For many applications, the beams’ initial transverse intensity profile has a single, on-axis peak such as a Gaussian, a ‘flat-top’ generated by a circular aperture, or some approximation thereof. Any additional peaks or non-monotonic decreases in the profile are often categorized as higher-order mode aberrations, and because they result in enhanced diffraction of the beam are usually undesirable [21]. Bessel and Airy beams are an exception. Bessel beams have cylindrically symmetric profiles with a central peak surrounded by concentric rings each possessing nearly the same amount of energy. During propagation the outer rings diffract inward fueling the center of the beam with energy and maintaining constant on-axis intensity. The inward diffraction also results in self-healing of the beam: when the beam’s propagation path is partially obstructed the profile can nearly reform itself down- range [22, 23, 24, 25, 26]. Airy beams, on the other hand, do not have cylindrical symmetry. The profile appears waffle-like with the intensity peaking in one corner and dropping with distance from the peak. While the Airy beam’s centroid follows a straight trajectory, the peak intensity propagates along a curved path. These curved paths approximately 21 follow the ballistic trajectory of projectiles and exhibit the same ‘stalling’ in motion that gravity produces in launched projectiles [26, 27, 28, 29]. As with the Bessel beam, the interference of diffracting beamlets in the Airy beam’s waffle-like pattern results in self-healing [26]. True Bessel and Airy beams are, however, only theoretical. Like plane waves, their field profiles would need to extend to infinity, requiring an unphysical, infinite amount of energy. In practice, the transverse profile of Bessel and Airy beams must be truncated at a finite radius or ‘sub-apertured,’ for instance with a Gaussian profile, changing how the beam propagates [23]. Nevertheless, even when apertured, the desirable features of Bessel and Airy beams survive over extended distances in vacuum. As a result, it has been suggested that Bessel beams [30] and Airy beams [26, 31] could be useful for remote sensing or directed energy applications. If the desirable properties of Bessel and Airy beams are to be harnessed for these applications, they must be robust to phase distortions resulting from atmospheric turbulence [32, 33, 34, 35]. In chapter 5 we investigate, through the use of propagation simulations, how turbulence modifies the propagation of Bessel and Airy beams. For both beams, we find the extent that each beam is modified by atmospheric turbulence depends on the transverse beam size. The transverse coherence length (Fried parameter), describing the transverse distance over which the phase fronts of the laser beam remain correlated, decreases with propagation distance. As a result, large diameter beams acquire transverse phase distortions that modify their propagation before they undergo standard diffraction. Small diameter beams, on the other hand, diffract before 22 their transverse coherence is disrupted. Furthermore we find that the nature of turbulence-induced beam spreading differs between Bessel and Airy beams and Gaussian beams: the simultaneous diffraction of many rings or beamlets make the Bessel and Airy beams resistant to spreading in turbulence. Despite the Bessel beams resistance to spreading, the power delivered to the initial aperture area decreases as the number of rings increases. This suggests that the most effective Bessel beam for delivering power is a Bessel beam with zero rings, essentially a Gaussian beam. 23 Chapter 2 : Reciprocity and Adaptive Optics 2.1: Overview Propagation of radiation through a medium is considered reciprocal if the relation between a point source and point receiver is invariant when their locations are interchanged [36]. Reciprocity is particularly significant when the relation between the point source and receiver cannot be explicitly derived, for example when considering propagation through a turbulent medium. This appears in applications such as directed energy or remote sensing, when it is desirable to focus a laser beam through atmospheric turbulence onto an uncooperative target. This is a difficult problem to solve because the state of the atmospheric channel is unknown, the target is located at an arbitrary location, and the radiation is spatially incoherent after scattering from the target. It has been demonstrated that adaptive optics can correct for the distortions caused by atmospheric turbulence in situations where the reflector is coherent or cooperative [37, 15, 38, 39, 16]. Adaptive optics relies on the reciprocity of the turbulent atmosphere [40], and for this reason there has been much work devoted to reciprocity in the field of atmospheric optics [15]. We begin our work on active optics in section 2.2 by using reciprocity to derive a relation between the field incident upon a reflector and the field detected by a point receiver. Then in section 2.3 we illustrate the derived relation in terms of a random walk in the complex plane. The analysis in section 2.3 suggests that the relation could 24 be useful in an adaptive optics system. After verifying that our numerical simulation agrees with the derived relation in section 2.4, we simulate an adaptive optics system which utilizes the relation to focus a laser beam through atmospheric turbulence onto a rough surface in section 2.5. While the relation is general to any reciprocal optical configuration, we discuss it here in the context of focusing a laser beam through atmospheric turbulence onto a target that provides a spatially incoherent return. Figure 2.1 displays the reciprocity preserving geometry considered throughout the chapter. The beam starts at the transmitter plane propagates through the turbulent channel, and undergoes scattering at the target plane. The scattered beam is collected in the receiver plane, where it is focused to the detector. The transmitter and receiver planes coincide, creating a monostatic channel. An aperture and phase modulator are located in the transmitter/receiver plane. A thin lens with focal length f is centered a distance f from the detector and transmitter/receiver planes. Figure 2.1: (a) Schematic of the propagation geometry. The transmitter and receiver plane coincide, creating a monostatic channel. The medium between the transmitter/receiver plane is arbitrary, but in this context we consider a turbulent atmosphere. A thin lens is located 25 between the receiver and detector planes, equidistant from each plane. (b) Schematic of the phase modulator located in the transmitter/receiver plane. The phase modulator consists of square elements of length a which apply piston type phase perturbations to the transmitted and received beams. The total dimension of the phase modulator is defined as b. 2.2: Derivation of the reciprocity relation As discussed in chapter 1, we express the transverse electric field, E , as a carrier wave modulating a slowly varying envelope, E : 1 2 ( , ) ( )exp[ ( )] c.c.E t E i kz t   x x where /k c . Since we are only interested in the electric field envelope at specific planes, we adopt the shorthand ˆ( )DE r , ( )OE r , ( )RE r , and ( )TE r for the complex electric field envelopes at the detector, transmitter, receiver, and target respectively, where rˆ , r , and r are the associated transverse coordinates. The spatial coordinates and fields are normalized by 1/2(2 / )f k and 1/20(2 / )OI c respectively, where OI is the peak intensity of the outgoing beam. The function ( ) r describes the transfer function of the aperture and phase modulator in the transmitter/receiver plane, which can be any complex, square integrable function. The initial envelope of the transmitted beam can be expressed as the field of a point source located at ˆ0r propagated from the detector plane to the transmitter plane and filtered by the function ( ) r : ˆ( ) ( )exp( 2 )OE i i    0r r r r . (2.1) We denote the point spread function, or Green’s function for propagation from the transmitter to the target as ( ; )TOG r r . In terms of the Green’s function, the electric field envelope at the target is then given by 26 ( ) ( ; ) ( )T TO OE G E d  r r r r r , (2.2) which upon making use of Eq. (2.1) provides ˆ( ) ( ; ) ( )exp( 2 )T TOE G i d     0r r r r r r r . (2.3) The target reflects the field, and in general both apertures the beam and imparts a phase shift, such that the reflected field has the initial profile ( ) ( )exp[ ( )]TE F i  r r r , where F is a real square-integrable function and  is real. Denoting ( ; )RTG r r as the Green’s function for propagation from the target to the receiver, we have ( )( ) ( ; ) ( ) ( ) iR RT TE G E F e d       r r r r r r r . (2.4) We stress that ( ) r can represent any phase shift imparted upon reflection. Back at the receiver plane, the field passes through the aperture and phase modulator and propagates from the front focal plane of the lens to the back focal plane, coinciding with the detector plane. The resulting envelope at the detector is then ˆ ˆ( ) ( )exp( 2 ) ( )D RE i i E d    r r r r r r . (2.5) Using Eq. (2.4) we can rewrite this as ( )ˆ ˆ( ) ( ) ( ) ( ; ) ( )exp( 2 )iD T RTE i E F e G i d d               r r r r r r r r r r r . (2.6) If the channel between the transceiver and target is reciprocal, the Green’s function for propagation from the transmitter plane to the target plane is equal to the Green’s function for propagation from the target plane to the receiver plane [40]: ( ; ) ( ; )TO RTG G r r r r . If we observe the field at the location ˆ0r we can write Eq. 27 (2.6) in terms of the envelope at the target plane described by Eq. (2.3), which gives the final result 2 ( )ˆ( ) ( ) ( ) iD TE E F e d      r 0r r r r . (2.7) Recall that the location ˆ 0 r is the location of the point source corresponding to the transmitted beam and can be determined by the tilt of the transmitted beam. Experimentally, this location can be found through the enhanced backscatter (EBS) phenomenon. An in-depth study on EBS is presented in chapter 3. Equation (2.7) relates a single quantity measured by the detector, ˆ( )DE 0r , to the field at the target plane ( )TE r , with information about the optical configuration contained only implicitly in ( )TE r . The remaining factor, ( )exp[ ( )]F i r r , captures the interaction with the target, and can take any form. Reflection from a mirror at normal incidence, for instance, is described by ( ) 1F  r and ( ) 0  r , while an idealized glint could be described by ( ) ( )g gF A    r r r , where gr is the origin of the glint on the target and gA its effective area. Of primary interest here is reflection from a rough surface, which can be described by a pointwise random phase ( ) ~ unif[0,2 ] r with ( ) 1F  r . At first glance, it might appear that application of such a phase would eliminate implicit information of the optical configuration contained in Eq. (2.7): with 1TE  , ˆ( )DE 0r averages to zero. However, an interpretation of Eq. (2.7) in terms of a random walk in the complex plane reveals its significance. 28 2.3: Random walk in the complex plane We conceptualize the continuous target surface as a grid of discrete scatterers. Each scatterer is separated by a distance l, occupies an area l 2 , and has a uniformly distributed random phase. This is equivalent to approximating Eq. (2.7) as a Riemann sum: 1 ˆ( ) n N i D n n E e   0r . (2.8) where 2 2| ( ) |n T nE l  r , 2arg[ ( )] ( )n T n nE   r r , nr is the location of the n th scatterer, N is the number of scatterers illuminated by the incident radiation, and we have mapped the 2D coordinate system to a 1D array for convenience. Equation (2.8) is a sum over random phasors, and can thus be interpreted as a random walk in the complex plane [41]. Typically the transverse coherence length of optical phase fluctuations [11] far exceeds the correlation length of a rough surface [42]. We thus ignore the phase contributed by the incident field in the following discussion such that unif[0,2 ]n  . Simulations presented below provide further justification for this. The magnitude of each step in the random walk is proportional to the intensity at a single location in the target plane, n , while the direction of each step is determined by n . The walk’s result is the field measured in the detector plane. For simplicity, we consider the uniform illumination of an incoherent reflector, corresponding to an isotropic, uniform random walk. This fixes all n to a single value, n  , such that the expected intensity in the detector plane is given by 2 20ˆ| ( ) |E DI E N r , where 29 the brackets denote the expected value. As a side note, the proportionality EI N can be understood as the result of an incoherent summation: , exp[ ( )] N n m n m N i    . (2.9) Because n is determined by the target, it cannot be modified, and the random walk will remain isotropic. However, by adjusting the outgoing beam’s phase, the intensity distribution on the target, i.e. the length of each step in the walk, can be modified. In particular, if the beam power, P N , is fixed, then const.N  , and EI P . This proportionality demonstrates that, on average, an intensity increase at the detector point 0rˆ translates to an increase in the target plane intensity. Figure 2.2 provides a simple example of this concept. One instance of an isotropic, random walk in the complex plane is displayed in Figure 2.2a with 64N  steps and 1  . The walk begins at the origin denoted by the ‘O’ and ends at the ‘X’. Each step in the random walk is denoted by a black arrow, while the overall result of the walk is denoted by the blue arrow. The expected length of a random walk with these properties is 8, and for this particular instance the length of the resulting walk is 5.7. On the other hand, Figure 2.2b displays a random walk with fewer steps, but each step is longer. Specifically there are 16N  steps each of length 4  . The quantity P N is the same for both walks; however the walk in Figure 2.2b has an expected length of 16. The specific walk in Figure 2.2b results in a length of 12.4. By redistributing the energy into less steps of greater magnitude, the expected length of the walk is increased. 30 Figure 2.2: (a) Isotropic random walk in the complex plane with 64 steps of length 1. (b) Isotropic random walk in the complex plane with 16 steps of length 4. 2.4: Numerical verification of reciprocity relation We continue by applying this analysis to focusing a laser beam through atmospheric turbulence onto a spatially incoherent target. Specifically, we simulate the optical configuration depicted in Figure 2.1a. A 2 / 1 mk    , collimated laser beam is initialized in the transmitter plane with a transverse Gaussian profile of spotsize 0w . Between the transmitter and the target, the beam evolves according to the steady-state paraxial wave equation with turbulence-induced refractive index fluctuations included through phase screens. We use the analytic approximation of the Hill spectrum, Eq. (1.22) , to characterize the distribution of turbulent scale sizes with an inner scale of 1 mm and outer scale of 1 m [1]. To simulate scattering from the rough surface at the target plane, we first apply a random phase uniformly distributed from zero to 2π to the electric field amplitude at each grid point [41]. We then apply a circular filter in the Fourier domain that 31 eliminates transverse wave numbers near the domain boundary. Propagation between the receiver and detector plane is modeled using a Fourier transform and appropriate scaling of the domain lengths based on the focal length of the lens. Figure 2.3: . Scatter plot of the magnitude (a) and phase (b) depicting the numerical equivalent of the relation described by Eq. (2.7). The vertical and horizontal axes correspond to the left and right hand sides of Eq. (2.7) respectively. Since the values are complex, there are two scatterplots to demonstrate that the simulation in good agreement with the theory. As a demonstration of Eq. (2.7) and to verify the numerical model, we simulated propagation through 1000 realizations of the optical configuration shown in Figure 2.1a. The turbulent channel was 2.5 km long, and composed of 10 phase screens each with a refractive index structure function, 2 nC , chosen randomly from a uniform distribution on the interval 16 13 2 3[10 ,10 ] m   . Choosing 2nC randomly for each realization provides verification over a wide range of turbulence-induced beam distortions in both the detector and target planes. The aperture function was set to a Gaussian, 2 2 0( ) exp( )r w  r , with 0 5 cmw  . Figure 2.3 displays scatter plots of the resulting left and right hand sides of Eq. (2.7) using the previously described normalization. Both the phase and amplitude predicted by the simulation are in 32 excellent agreement with Eq. (2.7). We note that the simulation includes the turbulence phase contribution neglected in the random walk discussion above. 2.5: Adaptively focusing through turbulence The random walk interpretation of Eq. (2.7) suggested a positive correlation between the focal and target intensities. Thus a feedback loop where the focal intensity informs the phase imparted to the transmitted beam, should, on average, result in an increased intensity on target. In the remainder, we describe simulations of an adaptive optics implementation that does just that. The adaptive optics system consists of a phase modulator (Figure 2.1b) in the transmitter/receiver plane controlled by the stochastic parallel gradient descent algorithm (SPGD), which we return to below. The phase modulator is comprised of square elements length a on a side in a 15 15 array. Each element can make piston phase changes. The elements are indexed by l and m corresponding to the x and y dimensions respectively, and are defined on the interval [ 7,7] . For example, the indices ( 7, 7)l m    , ( 0, 0)l m  , and ( 7, 7)l m  correspond to the bottom left, center, and top right elements respectively. We use the conventional stochastic parallel gradient descent (SPGD) algorithm and notation [37]. The SPGD algorithm is designed to iteratively increase a metric by applying random perturbations to all control parameters and calculating the gradient of the metric. Motivated by the preceding discussion, our metric is ˆ| ( ) |DJ E 0r . The phases applied by each element of the phase modulator serve as the control 33 parameters. In particular, the phase applied by element l, m at the end of iteration n is defined as ( 1) , n l mu  calculated in three steps: ( ) ( ) , , , n l m l m l mu u S    , ( ) ( ), , , n l m l m l mu u S    , and  ( 1) ( ) ( ) ( ), , ,n n n nl m l m l mu u S J J     . (2.10) where ,l mS is randomly chosen as +1 or -1 for each element,  is the magnitude of the perturbation, ( ) / nJ  corresponds to the metric value resulting from ( / ) ,l mu   , and  is the gain coefficient. For the simulation, we chose the adaptive optics parameters a = 2 cm, 10  , and 46 10   . The turbulent channel was 2.5 km in length with a 2 15 2 310 mnC   . This is considered mild optical turbulence; we elaborate on the motivation for these parameters and limitations of the adaptive optics system below. The transmitted beam had 0 10.6 cmw  clipped at diameter 15 30 cmb a  . The corresponding transfer function for the aperture and phase modulator during iteration n is given by  ( ) 2 2 ( )0( ) 1 ( 2) exp ( ) n nH r b r w iu       r r , (2.11) where H is defined as the Heaviside step function and ( ) ( )nu r is the phase applied by the phase modulator. Simulations were conducted for 300 statistically independent turbulent channels. Each simulation completed 100 iterations of the SPGD phase correction algorithm. Figure 2.4a and Figure 2.4b show the intensity profiles during iterations 1 and 100 respectively where the correction algorithm was most effective, i.e. producing the highest maximum intensity on the target. In this instance, the maximum intensity during the first iteration is 1.7 due to minor scintillation effects, but there does not exist a localized region of high intensity. After 100 iterations of the SPGD algorithm, 34 the metric increased from 304J  to 1,280J  . The beam during iteration 100, displayed in Figure 2.4b, exhibits a single localized region of high intensity, 9.5 times the transmitted beam’s maximum intensity with approximately 10% of the transmitted power contained in a radius of 1 cm about the peak intensity. Out of the 300 simulations there were only 9 instances where the correction algorithm did not increase the maximum intensity on the target. Figure 2.4c and Figure 2.4d display the intensity profiles during iterations 1 and 100 respectively for one such instance. In this simulation, the SPGD algorithm increased the metric from 77J  to 880J  , while the maximum intensity decreased from 2.6 to 2.1. Instead of a single, localized region of high intensity, there are multiple local intensity maxima. By chance, the incoherent sums for each of these maxima added constructively, and resulted in an increased metric without an increased intensity on the target. In practice, the combination of a moving target with a finite photodetector exposure time would effectively perform an average and may mitigate such an occurrence. 35 Figure 2.4: Top row: Initial (a) and final (b) intensity in the instance when SPGD was most effective. The metric increased from 304J  to 1,280J  bewtween (a) and (b). The intensity profile (a) before the SPGD algorithm is distorted by atmospheric turbulence. After 100 iterations of SPGD the intensity profile (b) displays a localized region of high intensity. Bottom row: Initial (c) and final (d) intensity profiles in the instance when SPGD was least effective. The metric increase from 77J  to 880J  between (c) and (d), however there exist multiple local intensity maxima in (d) instead of a single high intensity region. Statistics over all 300 simulations are displayed in Figure 2.5. Figure 2.5a and Figure 2.5b display histograms of the maximum intensity and power in the bucket (PIB) on the target. The PIB is calculated as fractional power contained within a 1 cm radius bucket centered on the location of maximum intensity. Statistics obtained from intensity profiles during iterations 1 and 100 are plotted as the red and blue bars respectively. Due to minor scintillation, the average maximum intensity during 36 iteration 1 is 1.9 with a standard deviation of .36. After 100 iterations of SPGD the average maximum intensity is 3.5 with a standard deviation of 1.2. The same trend is seen in the PIB. The average PIB during iteration 1 is .027 with a standard deviation of .005 and the average PIB during iteration 100 is .045 with a standard deviation of .015. Figure 2.5: (a) Histogram of maximum intensity over 300 independent channels comparing the intensity profiles at iteration 1 (red) and iteration 100 (blue) of the SPGD algorithm. The average maximum intensity during iteration 1 is 1.9 with a standard deviation of .36. The average maximum intensity during iteration 100 is 3.5 with a standard deviation of 1.2. (b) Histogram of PIB contained within a 1 cm radius bucket centered on the location of maximum intensity. The average PIB during iteration 1 is .027 with a standard deviation of .005. The average PIB during iteration 100 is .045 with a standard deviation of.015. 2.6: Discussion on reciprocity and adaptive optics In choosing the parameters for adaptive focusing based on Eq. (2.7), several conditions and limitations must be considered. First, there cannot be substantial intensity scintillation. Since our adaptive optics implementation only applies phase corrections, it cannot effectively correct for strong scintillation [43, 44, 45]. For weak turbulence, the scintillation index equals the Rytov variance, 2 2 7 6 11 61.23R nC k z  , providing the condition 2 .2R  . For our turbulent channel, the Rytov variance was 37 2 .18R  . Recent research in combined amplitude and phase corrections may relax this condition [45]. Second, the dimension of the phase modulator elements must be less than the transverse coherence length, 0a  , where 2 2 3 5 0 1.67( )nk C z  . Ideally a fast algorithm for correcting the wavefront would allow arbitrarily small phase modulator elements. However, with iterative algorithms, decreasing the element dimension can require decreasing the perturbation magnitude,  , which increases the convergence time. Third, iterative algorithms find local optima, which are not necessarily the ideal global optimum [37]. More sophisticated algorithms may offer a solution, but are a fundamental adaptive optics problem outside the topic of this dissertation. Finally, the transit time of the light to the target and back must be shorter than the time scale over which the turbulence changes. Turbulent fluctuations typically become uncorrelated over a time scale of milliseconds, which does not present a problem for typical propagation distances. However this time scale depends on a variety of parameters such as wind, temperature, and humidity; and in extreme instances the turbulent atmosphere may not possess the property of reciprocity. Working within these limits, we have demonstrated the derived relation can facilitate the formation of a localized region of high intensity on a remote, spatially incoherent target. The derived relation is quite general, relying entirely on reciprocity, making it suitable for applications requiring propagation through random media. 38 Chapter 3 : Enhanced Backscatter 3.1: Background When an optical beam propagates through turbulence to a target, reflects off the target, and propagates back through the same turbulence, it undergoes a local recovery of spatial coherence and intensity enhancement referred to as enhanced backscatter (EBS) [1, 46, 47, 48, 49, 50, 51, 52, 53, 54]. EBS is observed as a retroreflected, localized intensity spot in the detector or receiver plane when many independent instances of turbulence are averaged over. The EBS location coincides with the location 0ˆr in the previous chapter, providing an easy method of accurately finding the optical axis of the system. EBS was extensively studied after the first theoretical predictions of the phenomenon in 1972 by Belenkii and Mironov [55]. Numerical investigation followed shortly after using the phase screen model to represent turbulence induced phase distortions [52]. EBS was later observed in lab experiments [53] and field experiments [54] in 1977 and 1983 respectively. Further studies have demonstrated the potential of EBS for precision pointing and tracking [56]. However, there have been fewer recent investigations, possibly due to minimal applications and detection difficulty, but major advances in both computing and imaging technologies warrant further studies of EBS. The objective of this chapter is fivefold: First to experimentally measure EBS from retroreflectors and diffuse reflectors using modern imaging systems. Second, to 39 validate propagation simulations of EBS with experimental results. Third, to demonstrate a scaling between the experiments and DE-scales, justifying use of the simulations for examining EBS in this regime. Fourth, to verify theoretical predictions on the influence of turbulence strength on EBS. And finally, to present a novel method for detecting EBS, which greatly reduces the time required for observation. We use a similar geometry to Figure 2.1 which is displayed in Figure 3.1. As in chapter 2, we define the transmitter, target, and receiver planes as the planes where the beam starts, reflects, and ends respectively and a double pass refers to the beam propagating from the transmitter plane to target plane and back to the receiver plane. A bistatic channel refers to a channel where the turbulence on the path to the target is uncorrelated with the turbulence of the return path. Bistatic channels arise when the reflected beam path does not overlap with the initial beam path or when the turbulence correlation time is shorter than the double pass transit time. A monostatic channel refers to a channel where the turbulence is identical on the path to the target and the return path. Monostatic channels arise when a beam is reflected at a small angle back toward the transmitter and the turbulence correlation time is longer than the double pass transit time. EBS can be understood through a ray model of a collimated beam propagating through a turbulent channel and reflecting from a diffuse reflector. In the ray model, a collimated laser beam consists of many individual rays starting with zero transverse wavenumber in the transmitter plane. As the rays propagate through the atmosphere, the turbulent fluctuations in the refractive index change each ray’s transverse 40 wavenumber. When the rays reach the target plane they are scattered back toward the transmitter plane along different paths than those taken to the target. A subset of the rays will end with zero transverse wavenumber in the transmitter plane and within the aperture of the outgoing beam. If the channel is monostatic, because these rays end up within the aperture of the outgoing beam, there must exist a ray that started at the same location and took the same path through the turbulence but in the opposite direction. These rays are referred to as reciprocal rays and always come in pairs. Reciprocal rays pass through the same turbulence cells but in reverse order. Since their path length is the same, the two rays forming a reciprocal ray pair are always in phase at the receiver plane [47]. The reciprocal rays will be normally incident on a lens placed in the receiver plane, and refracted to the focal point. While a reciprocal ray pair is always in phase, each pair is generally out of phase with all other reciprocal pairs: each pair has a different path length through the turbulence. This leads to high average intensity and a high variance in intensity at the focal point of the lens, known as enhanced backscatter. The intensity at the focal point is in general predicted to be increased by a factor of two [47, 57]. The ray model, displayed in Figure 3.1, gives a concise, visual overview of EBS, yet provides a remarkable level of insight into the phenomenon [47, 57, 56]. For a comprehensive wave theory derivation of EBS see [48] and [50]. 41 Figure 3.1: Diagram of a pair of reciprocal rays taking the same path in opposite directions. EBS is a consequence of the reciprocity of atmospheric turbulence [40]. In a monostatic turbulent channel, the point in the detector plane that possesses the properties described in chapter 2 coincides with the location of the EBS peak. Thus observation of the EBS peak reveals a point at which, in theory, instantaneous measurements of the beam intensity on target can be made [38]. Furthermore, observation of EBS ensures that the atmospheric channel is reciprocal, which is not true when long distances or abnormally high wind velocities are considered. This is our motivation to investigate EBS. In this chapter we examine EBS on two length scales. The smaller scale corresponds to lab experiments using retroreflector and diffuse reflector targets in which fan assisted hot plates generate strong turbulence over a short distance. We show that the standard simulation model for propagation through turbulence, a steady-state paraxial wave equation with optical turbulence modeled with phase screens [1, 2], largely agrees with experimental observations of the EBS phenomenon. We then use the simulation to consider a second, larger length scale 42 relevant to DE applications. At these length scales, we find that EBS detection requires weak and strong levels of turbulence for retroreflectors and diffuse reflectors, respectively. We also find that although the EBS phenomenon is robust at this larger scale, the standard method of averaging over independent instances of turbulence is not practical. Performing a temporal average over turbulence causes two problems: a delay in detecting the EBS signal, and when the EBS signal is detected, it does not represent the current state of the channel. We instead introduce a novel method of EBS detection called the tilt-shift method (TSM), in which a spatial average replaces the standard temporal average. We begin in section 2 by describing the experimental arrangement for measuring the intensity profile of a laser beam propagated through lab-scale, hotplate generated turbulence. We then verify the propagation simulation reproduces the experimental results. In section 3, we scale the simulation to parameters relevant to DE applications and present results comparing the EBS signal from retroreflectors and diffuse surfaces. Section 4 describes the TSM which enables the detection of EBS in frozen turbulence. Section 5 ends the chapter with a summary of our results. 3.2: Lab-scale turbulence Experimental arrangement A schematic of the experimental arrangement is shown in Figure 3.2 for both the monostatic and bistatic cases. An initially collimated beam with wavelength 635 nm  and spot size 0 1 cmw  originates from a Thor Labs fiber laser. The beam propagates through hotplate generated turbulence, reflects off the target then 43 propagates back through the turbulence once more before being focused onto the CCD camera. The same optical components are present in the bistatic and monostatic channels so that any loss or deformations to the beam from the components is the same for each channel. Figure 3.2: Schematic of the experimental arrangement for the (a) monostatic and (b) bistatic channels. Turbulence is generated by two 61 cm hotplates for each channel. The turbulent channel is 3 m long and fueled by two 60 cm 7.5 cm hotplates that can reach up to 200 o C. Through thermal conduction the hotplates generate a turbulent boundary layer. The addition of a duct fan provides a jet stream down the axis of the channel that facilitates energy cascade and sustains the turbulent flow. The hotplates are wide enough to allow a slant pass across the channel for the transmitted beam in the bistatic arrangement. 44 The experiment was conducted with two types of target: a retro and a diffuse reflector. A 2 inch diameter dielectric cornercube was used as the retroreflector, and a simple white poster board for the diffuse reflector. The reflected beam passed through a 50/50 beam splitter and was focused with a 50 cm focal length, 2 inch lens (5.08 cm) onto a CCD camera with 1024 1024 square pixels 5.5 µm on a side placed in the focal plane. The exposure time of the camera was set to 30 μs and 15 ms for the retroreflector and diffuse reflector respectively. The intensity was averaged over 100 frames collected at 8 frames per second. Simulation description As described in chapter 1, we use the steady-state paraxial wave equation with optical turbulence modeled with phasescreens to solve for the transverse electric field of a laser beam propagating through the turbulent channel.. The transverse simulation domain is 204.8 mm  204.8 mm with a resolution of .1 mm for the retroreflector, and 51.2 mm  51.2 mm with a resolution of .025 mm for the diffuse reflector. The turbulence is characterized by an analytic approximation to the Hill spectrum, Eq.(1.21). Following the method described in [58], we experimentally measured the scintillation index and mean squared angle of arrival to be 2 4(6.3 .2) 10I    and 2 8 2(1.4 .1) 10 rad    respectively. Using these values along with the equations [59]: 2 2 1 3 03.28 nC Ll  , (3.1) 2 2 3 7 3 012.8I nC L l  , (3.2) 45 we calculate the inner scale and effective 2 nC associated with our hotplate generated turbulence to be 0 8.7 .5 cml   , and 2 10 2 3 , (6.2 .7) 10 mn effC     respectively. The outer scale is chosen to be 0 60 cmL  corresponding to the length of one hotplate. Fluctuations larger than the simulation domain are not accurately sampled. However, we are primarily concerned with small scale fluctuations that lead to scintillation so this does not present an issue. We note that this could be resolved through the addition of sub-harmonics in the process of generating phasescreens [60]. A double pass in a bistatic turbulent channel uses eight phasescreens: four for propagation from the transmitter plane to the target plane, and four for propagation from the target plane to the receiver plane. A double pass in a monostatic turbulent channel uses four phasescreens: the same phasescreens are used for propagation in both directions. The m th phasescreen is placed at a distance  1 2m L n where n is the number of phasecreens and L is the distance from the transmitter to the target. We model the propagation from the receiver plane to the detector plane, the focal plane of the lens, by taking the Fourier transform of the electric field in the receiver plane and scaling the domain lengths to correspond to focusing by a 50 cm focal length lens. All plotted quantities are considered in the detector plane. Intensity averages are performed over 1000 runs through statistically independent turbulent channels. Retroreflector First we examine the experimental and simulation results for the cornercube target. The cornercube inverts and reverts the wavefront corresponding to a 180 degree rotation and a reversal in the direction of propagation. This results in many reciprocal 46 ray pairs leading to a distinctive EBS intensity peak in the detector plane. Figure 3.3a and Figure 3.3b display the average intensity profiles for the monostatic and bistatic channels respectively obtained using the experimental arrangements described in section 2. Figure 3.3c compares the y = 0, x-slice for both channels. Intensity is normalized according to the equation ( , ) ( , ) (0,0)norm biI x y I x y I where (0,0)biI is the intensity at the center of the bistatic beam. The average intensity in the monostatic channel exhibits the sharp, central peak reaching a value of ~2 indicative of an EBS enhancement. Other than the EBS enhancement, the monostatic and bistatic channels produce nearly identical Gaussian intensity profiles. Figure 3.3: Average monostatic (first column), bistatic (second column) and x-slice (third column) intensity profiles from a retroreflector target obtained from the experiment (first row) and simulation (second row). The experimental and numerical results are in good agreement. The localized region of high average intensity is observed in the monostatic channel and is increased by approximately a factor of two in comparison to the bistatic channel. 47 In the simulation, the retroreflection at the target plane is applied by rotating the propagation domain 180 degrees and propagating the beam in the opposite direction. Figure 3.3d, Figure 3.3e, and Figure 3.3f display the average intensity profiles obtained from the simulation for the monostatic channel, the bistatic channel, and the y = 0, x-slice for both channels respectively. The turbulent beam spreading in the simulation matches the experiment, indicating that the strengths of turbulence are approximately equal. The primary difference in the intensity profiles is the sharper contrast of the EBS peak to the background predicted by the simulation. This could be due to small errors in alignment or aberrations from the optical components in the experimental arrangement, as the simulation assumes perfect alignment and flawless optical components. Regardless of this small difference, the 2-dimensional phasescreen simulation accurately models the 3-dimensional turbulence of the experiment and reproduces the EBS observed from a retroreflector. Diffuse reflector We now examine the experimental and simulation results for a diffuse reflector. In the experiment a large white poster board provides the diffuse reflection. Figure 3.4a, Figure 3.4b, and Figure 3.4c display the observed average intensity profiles for the monostatic channel, the bistatic channel, and the y = 0, x-slice for both channels respectively. As before, the intensity is normalized to the value at the center of the bistatic beam. The reflection from the rough surface creates a broad speckle pattern in the intensity profile for both the monostatic and bistatic channel. Again the monostatic channel produces an EBS peak with twice the bistatic intensity, but the peak is much wider than that resulting from the cornercube target. 48 Figure 3.4: Average monostatic (first column), bistatic (second column) and x-slice (third column) intensity profiles from a diffuse reflector target obtained from the experiment (first row) and simulation (second row). In the experimental arrangement the target is a flat, white poster board. The EBS peak is approximately twice as wide in the experiment compared to the simulation results. This is likely due to the material properties associated with the target [61]. A diffuse reflection is difficult to model precisely due to both the random and fractal nature of the surface when viewed on the scale of optical wavelengths. The surface scatters the beam to large angles contributing large transverse wavenumbers, | |k , to the electric field, which, in turn, cause significant beam spreading. Thus it becomes numerically intractable to make the Fourier and real space domains large enough to both resolve the large | |k values and follow the beam spreading. However, we are able to ignore | |k values larger than our Fourier domain for two reasons. First, the resulting reduction in beam power does not significantly affect the average intensity ratio of beams that propagate through different channels but are 49 reflected from the same surface. Second, the small scale phase distortions occur at random uncorrelated points and their associated | |k ’s leave the propagation domain of interest upon diffraction. We simulate a diffuse reflection in the same manner as in chapter 2. We first apply a random phase uniformly distributed from zero to 2 to the electric field envelope at each grid point [62]. The random phases are pairwise uncorrelated. We then apply a filter in the Fourier domain that truncates 5| | 1.1 10 rad m  k to ensure that all directions are limited by the same | |k . This step is required because the simulation domain is square. We note that if the reflection model creates a | |k that is larger than our Fourier domain, it will be aliased to another location within our Fourier domain—a property of discrete Fourier transforms. Normally this would corrupt a Fourier transform, but since the phase distortions from the diffuse reflection are random and uncorrelated, the aliasing results in a slightly different but equivalent random phase. Comparison of low-pass filters with different cutoff frequencies reveals that the cutoff frequency has no discernable impact on the EBS peak. Figure 3.4d, Figure 3.4e, and Figure 3.4f display the average intensity profiles obtained from the simulation for the monostatic channel, the bistatic channel, and the y = 0, x-slice for both channels respectively. The monostatic and bistatic intensity profiles both display a speckle pattern characteristic of a diffuse reflection as displayed in Figure 3.4a and Figure 3.4b. The EBS intensity enhancement factor of 2 is observed in the monostatic channel simulation. There is, however, a disagreement in the width of the simulated and measured EBS peaks. In particular, the measured EBS peak is approximately twice as wide as the peak predicted by the simulations. 50 This is most likely due to differences in the actual reflections from the white poster board target and those resulting from our diffuse reflector model. It has been shown that the roughness and material properties of the scattering surface impact the width of the EBS peak [61]. Since we are investigating EBS in the context of DE where we are generally unable to know the target’s material properties, we are not concerned with the difference in EBS peak width. 3.3: Directed energy scale Transition to DE length scales Using the Rytov variance and the Fried parameter [1]: 2 2 7/6 11 6 01.23R nC k L  , (3.3) 2 2 3 5 0 01.67( )nr C k L  , (3.4) as measures of the turbulence strength and transverse coherence length at DE scales, respectively, we can scale the simulation parameters for larger beams propagating over longer turbulent channels. This allows us to investigate EBS on length scales relevant to common DE parameters while avoiding the procedural difficulties of conducting experiments at this scale. First we calculate the Rytov variance and transverse coherence length for the lab- scale turbulence. Using the Rytov variance purely as a measure of turbulence strength, our parameters give the value 2 .83R  . Due to the large inner scale and short propagation distance of the lab generated turbulence, we use a different form for the transverse coherence length which depends on the size of the inner scale [1]:   1 2 2 2 1 3 0 0 01.64 nC k Ll   . (3.5) 51 We now demonstrate that a lab-scale and DE-scale simulation are nearly identical up to a scale length when the Rytov variance and the ratio of the transverse coherence length to initial beam radius are fixed. We define the beam size ratio as DE labw w  where DEw and labw are the 1 e radii of the electric field for the DE scale and lab-scale respectively. Noting the relation between the Rytov variance and Fried parameter in Eq.(3.3) and Eq.(3.4): 5 6 2 2 0 0 ~R L r k        , (3.6) we obtain the following expressions for the equivalent length and 2 nC of the DE-scale channel: 2 6 5 2 0 0.3( ) ( )RL k  , (3.7) 2 2 6 5 11 3 3 0 07.4( ) ( )n RC k     . (3.8) Using the values of the experiment 2 .83R  , 0 1.2 mm  , 1 cmlabw  , and 20 cmDEw  , we find 1.4 kmL and 2 15 2 38 10 mnC   for the equivalent DE- scale channel. The inner scale for the DE channel is chosen to be 1 mm, a common value for atmospheric turbulence. 52 Figure 3.5: Comparison of lab and DE simulations for a retroreflector (a) and diffuse reflector (b). With the appropriate scaling, we are able to find equivalent large-scale parameters which correspond to the lab-scale experiment. Figure 3.5a and Figure 3.5b compare the average monostatic intensity profiles of the lab-scale and equivalent DE-scale simulations for a retroreflector and diffuse reflector respectively. The x-axis of the large-scale intensity profile has been dilated by the factor 20  . Although the inner scales differ by approximately two orders of magnitude, the average monostatic intensity profiles nearly match. It is important to note that at the lab scale, both the beam width and final transverse coherence length are smaller than the inner scale. Scintillation develops when small and large scale phase fluctuations develop into amplitude fluctuations [63]. The small scale phase fluctuations arise from turbulent cells smaller than the transverse coherence length or beam width, whichever is smaller. However, the large inner scale of the laboratory generated turbulence greatly diminishes the amplitude of the small scale fluctuations. As a result, the scintillation observed in the focal plane results primarily from the lens-like, refractive phase perturbations imparted by the large scale turbulent fluctuations [63]. This is in contrast to the DE scale where both the large and small scale phase fluctuations contribute to scintillation. In spite of this difference between the lab and DE scale, Figure 3.5 suggests that the EBS signal, in 53 the far field, is relatively insensitive to the exact scale sizes of the initial phase perturbations. Directed energy parameters Having compared the simulation data with experimental data at the lab scale, we now turn to parameters used in DE applications. We simulate the propagation of an initially collimated Gaussian beam with wavelength 1 μm  and spot size 0 10 cmw  . The transverse simulation domain is 2.048 m by 2.048 m with a resolution of 1 mm. The distance from the transmitter plane to the target plane is 1 km. The turbulence is characterized by the Hill spectrum with 0 1 mml  and 0 1 mL  , which are common parameters for the inner and outer scale of atmospheric turbulence [1]. We consider 2 nC values varying from 16 2 35 10 m  to 13 2 35 10 m  corresponding to weak and strong atmospheric turbulence respectively. A double pass in a bistatic channel uses eight phasescreens and a double pass in a monostatic channel uses four phasescreens. All other parameters remain unchanged. To quantify the intensity enhancement, it is convenient to define the enhancement factor, F, as m b I F I  , (3.9) where mI is the average intensity at the focal point in the detector plane of a monostatic channel and bI is the average intensity at the focal point in the detector plane of a bistatic channel. As previously explained, the expected value of the enhancement factor is two under conditions of sufficiently strong turbulence. 54 Retroreflector Figure 3.6a displays the value of the enhancement factor averaged over 400 independent turbulent channels as a function of 2 nC and the corresponding Rytov variance for a corner cube target. Figure 3.6b displays the enhancement factor values over the first 200 runs at 6 different values of 2 nC . Figure 3.6 illustrates two major points: EBS requires averaging over many double passes, and strong turbulence is beneficial for EBS detection whereas it is usually a hindrance to long distance propagation. Figure 3.6: Enhancement factor as a function of 2 nC (a) and number of independent runs averaged (b) for a retroreflector target. The trend of enhancement factor versus 2 nC is in good agreement with early theoretical predictions. Plot (b) demonstrates that many runs must be averaged over to observe the enhancement factor of 2. 55 At 2 16 2 31 10 mnC    the turbulence is too weak for the EBS phenomenon to be easily observed and the enhancement factor is approximately 1. As 2 nC increases to 15 2/35 10 m  ( 2 .17R  ) the enhancement factor increases to 2. With a cornercube target, reciprocal paths exist even in the absence of turbulence due to the wavefront inversion in the target plane. Therefore the turbulence is not required to laterally shift rays to form reciprocal pairs and the turbulence only needs to be so strong as to produce minimal scattering effects. The relatively low levels of turbulence required for backscattered enhancement from a cornercube makes these targets suitable for experimental investigations of EBS as found in section 2. Theoretical predictions of the enhancement factor dependence on turbulence strength divide the turbulence strength into three regimes: weak, strong, and saturated turbulence. EBS is first observed in the weak fluctuation regime, where the enhancement factor transitions from 1 to 2. In the strong fluctuation regime, the enhancement factor temporarily reaches a maximum of approximately 3. As the strength of turbulence increases into the saturated regime, the enhancement factor converges to 2 [50, 51]. Figure 3.6a exhibits the same general behavior. The fluctuations in the enhancement factor are due to averages at each 2 nC value being taken over 400 statistically independent simulations. The error in the enhancement factor, F , is estimated using propagation of errors where the error for the average monostatic and bistatic intensities are determined by their standard errors respectively. Specifically, we use 1/2 1 2 2 2 1/2( ) [ ( / ) ]b m b m bF I N I I      , where 2i is the variance for each channel. Strictly speaking, the result has not fully converged. However it is not 56 difficult to identify the weak, strong, and saturated regimes as 2 .07R  , 2.07 1R  , and 2 1R  respectively. Figure 3.6b confirms that the data in Figure 3.6a was sufficiently averaged. It also illustrates the necessity to average over many instances of turbulence to obtain an accurate enhancement factor. In particular, the enhancement factor for 2 145 10nC   temporarily drops below 1 even though the long term average is approximately 2. We note that these curves correspond to one specific set of statistically independent simulations, and a different set of runs with the same 2 nC will average in a different manner but produce a similar result. Diffuse reflector We repeat this investigation for a diffuse reflector with the same parameters. Unlike a retroreflector, a diffuse reflector in the absence of turbulence does not produce reciprocal ray pairs. Although some rays get reflected directly backward, these are not considered reciprocal ray pairs because each ray essentially pairs with itself. Therefore, any reciprocal ray path must be the result of turbulence induced scattering, and, as a result, a diffuse reflector requires stronger turbulence than a cornercube for EBS observation. 57 Figure 3.7: Enhancement factor as a function of 2 nC (a) and number of independent runs averaged (b) for a diffuse reflector target. The trend of enhancement factor versus 2 nC is in good agreement with early theoretical predictions. Plot (b) demonstrates that many runs must be averaged over to observe the enhancement factor of 2. This is reflected in Figure 3.7a which displays the value of the enhancement factor averaged over 400 independent turbulent channels as a function of 2 nC and the corresponding Rytov variance. The enhancement factor first reaches a value of 2 at 2 14 2 3 24.2 10 m , 1.4n RC      which is 20 times greater than the minimum value of 2 nC required for the same enhancement from a retroreflector target. Correspondingly, the weak, strong, and saturated regimes are shifted to 2 1.4R  , 21.4 10R  , and 2 10R  respectively. 58 Figure 3.7b displays the approach to convergence for the enhancement factor as a function of the number of independent runs averaged for 6 values of 2 nC . This figure illustrates an issue not encountered with a retroreflector. As indicated by Figure 3.7a, turbulence characterized by 2 15 2 31 10 mnC    is not within the EBS regime for a diffuse reflector. However Figure 3.7b shows that the turbulence is strong enough to produce a transient enhancement factor, approximately 2, when too few runs are averaged over. This transient enhancement factor is the result of scintillation and not EBS. Atmospheric turbulence typically becomes uncorrelated over a time scale on the order of milliseconds [64] due to atmospheric conditions such as wind, temperature, and humidity. This means that intensity profiles should be sampled at a rate in the hundreds of Hz and it would take too long in the context of DE applications to distinguish an intensity enhancement due to EBS and an intensity enhancement due to scintillation. 3.4: Tilt-Shift method The standard method of detecting EBS requires averaging the intensity in the detector or receiver plane over many double passes through statistically independent turbulence. This puts a lower limit on the time delay between successive passes to ensure the turbulence has become uncorrelated. Furthermore, once the EBS signal has been obtained, it is only representative of the turbulence’s average state not the current state. Here we introduce a novel method, which we refer to as the tilt shift method (TSM), that can detect EBS from a single instance of turbulence. As a result 59 the method eliminates the time delay, and provides an EBS signal representative of the current state of the turbulence. The method expands on the principle of reciprocity in the EBS phenomenon: the same distortions are applied to the beam on return path as the outgoing path. This implies we can modify the transmitted beam, as long as we make the same modification to the reflected beam. In the TSM, we apply a small angle tilt to the transmitted beam before it enters the turbulent channel. Instead of re-tilting the reflected beam after it exits the channel, we shift the coordinates in the detector plane according to the initial tilt. This is equivalent to re-tilting the beam because the detector plane is the focal plane of the lens. By shifting the coordinates, the reciprocal rays refract to the focal point while the majority of the reflected beam is diverted away from the focal point. We note that the coordinate shift can be accomplished by a processor and no extra hardware is required. Additional beams, each with a unique tilt, are propagated through the same channel, and given the appropriate shift in the detector plane. Each beam propagates along a different path and contributes a different set of reciprocal rays at the focal point. The EBS can then be observed in the intensity profile averaged over small angle tilts through the same turbulent channel. Said differently, the TSM replaces an average over paths through statistically independent channels for an average over different paths through the same channel. We note that the reflected beam path must overlap with the transmitted beam path, placing an upper limit on the tilt angle. While the tilt shift method works for retroreflector and mirror targets, we choose to focus on a diffuse reflector because the wide angle of reflection increases 60 the upper limit on the tilt angle. To simulate a tilt, we apply a complex phase to the initial field envelope of the i th beam, ,, ,( ,0) ( ,0) tilt ii i iE E e    k r r r where the magnitude of the tilt wavenumber, , 0| |tilt i tilt tiltk k  k , and absolute tilt angle, tilt , are equal for each beam. A unique tilt is defined through the angle i where , , cosx tilt i tilt ik k  and , , siny tilt i tiltk k  . The coordinates in the detector plane are shifted by , 0tilt i f kk where f is the focal length of the lens. Figure 3.8: Average monostatic (a) bistatic (b) and x-slice (c) intensity profiles using the tilt shift method in static turbulence. The region of localized high intensity is observed from a single instance of turbulence. This allows observation of EBS without performing a time average. Figure 3.8 shows the monostatic and bistatic intensity profiles from the TSM using 16 equally spaced tilts ( 1 / 8i i     ) of magnitude 52 tilt rad  . The DE parameters described in section 3 are used along with 2 14 2 35 10 mnC    . Each beam propagates through the same frozen turbulent channel. We observe EBS in the monostatic channel from a single instance of turbulence. 61 3.5: Discussion on enhanced backscatter We have observed EBS from both retroreflector and diffuse reflector targets. An effective value of 2 nC and inner scale were determined for the lab turbulent channel. These values were then used to verify that the phasescreen simulation agrees with the experimental measurements of EBS. We showed that the simulation can be scaled to length scales associated with DE applications while preserving the characteristics of EBS. Using the simulation, we investigated the dependence of turbulence strength on the enhancement factor. These results agreed with theoretical predictions. This investigation revealed that a diffuse surface requires much stronger turbulence than a retroreflector for observation of EBS. However, EBS is still a robust phenomenon for DE applications due to the long distances and high levels of turbulence commonly encountered. It also revealed that a rapid method of detecting EBS is required if EBS is to be used in DE systems. The TSM fulfills this requirement because it allows detection of EBS in frozen turbulence, thus eliminating the need to wait for the turbulence to change. 62 Chapter 4 : Beam Combining 4.1: Overview of beam combining Coherent combining has proven to be effective in situations with very low-power lasers and weak turbulence [19]. However, these conditions are not applicable for DE systems which require high-power lasers and must be effective in conditions of moderate and strong atmospheric turbulence. High-power lasers have significantly broader linewidths due to stimulated Brillouin scattering, Doppler shift, self-phase modulation, and Raman broadening in the gain medium. Stimulated Brillouin scattering (SBS) is the main contributor to the large linewidths in high-power fiber lasers. The SBS instability results in a threshold power level for fiber and slab lasers. To increase the threshold power level, the linewidth is intentionally broadened. The full width half max (FWHM) of the linewidth for high-power lasers (multi- kW) is typically 10GHz  and the coherence time is 1 0.1nsecCt   (for all single peaked lineshapes) [65]. Phase locking of lasers requires measuring the output phase information and applying the corrective phase to the individual lasers. To be effective, the time scale required for this process should be less than the laser coherence time. The coherence time is the characteristic time over which the phase and intensity randomly vary and is due solely to the finite spectral linewidth in the gain medium and the statistical nature of the emission from the atoms or molecules [4]. Stated simply, if the process of measuring the output phase and applying the corrective phase takes longer than the coherence time, then the output phase has 63 changed before the corresponding corrective phase can be applied. To the best of our knowledge, there do not currently exist instruments that can operate at rates comparable to the FWHM of the linewidth for high-power fiber lasers. However, in the following work we assume that there exist instruments that can control the phase over time scales shorter than 1 0.1nsecCt   . Furthermore, we assume that these instruments operate with zero error. As previously stated, the power spectral linewidth associated with the gain medium results in temporal fluctuations of the intensity and phase at the output. The nature of the temporal fluctuations depends on the shape of the power spectral linewidth. Fluctuations associated with a Lorentzian power spectral density are characteristically different from the fluctuations associated with a Gaussian power spectral density due to the extended wings of the Lorentzian spectrum. A Lorentzian spectrum can result from collisional broadening, while a Gaussian spectrum can result from Doppler broadening or local random atom environments in a solid. Regardless of the lineshape, the standard deviation of random fluctuations in the output intensity is approximately equal to the average intensity. In addition, the output phase of the field randomly fluctuates over 2 radians on time scales down to the coherence time. Experimental observations on the lineshapes of high-power fiber lasers have revealed lineshapes that are approximately Lorentzian. In this chapter we discuss the physical processes associated with the combining and propagation of high-power laser beams. We compare the energy delivered to a target for the case of coherently combined and incoherently combined laser beams. Through simulation, we demonstrate that in vacuum and weak turbulence, the 64 effectiveness of coherent combining is limited by the phase fluctuations occurring due to the broad linewidth of high-power lasers (Figures 7 and 8). Furthermore, the advantages of coherent combining begin to diminish in moderate turbulence, as the turbulence induced phase distortions become dominant (Figures 10 and 11). In strong turbulence and multi-km propagation ranges, we find negligible differences in the energy on the target between coherently and incoherently combined laser beams (Figures 12 and 13). The incoherently combined architecture is far simpler to implement and both beam combining architectures have the capability of employing adaptive optics to extend the range. We begin in Section 2 by describing the model for the beam director, large linewidth source, and coherent combining instruments. We present our propagation simulation results in Section 3 for vacuum, moderate turbulence, and strong turbulence. A discussion of our findings concludes this chapter in Section 4. 4.2: Beam Director Model Geometry of the beam director We model the beam director as a set of nine square beams which we refer to as tiles. The side length of each tile is a and the distance between the centers of adjacent tiles is b where b a . For example, a beam director with the parameters a=b would have no filling factor and appear as a single tile with side length 3a and centered at  0, 0x y  . The tile arrangement in the xy plane is shown in Figure 4.1a. Tiles are indexed by the indices l and m corresponding to the x and y dimensions respectively. Indices l and m are integers defined on the set  1,0,1 . For example, the indices 65  1, 1l m    ,  0, 0l m  , and  1, 1l m  correspond to the bottom left, center, and top right tiles respectively. Figure 4.1: Beam director geometry. We define the side length of each tile as a and the distance between the centers of adjacent tiles as b where b a . Each tile is indexed by the indices l and m defined on the set  1,0,1 . Figure 4.1b shows a side view of the beam director, i.e. yz plane. Each tile is positioned with tip-tilt correction so that the center of the tile propagates to the point  0, 0,x y z L   in vacuum where L is the distance from transmitter to the target and L b . The tilt applied to a tile with index m is defined as y mb L  . Similarly the tip applied to a tile with index l is defined as x lb L  . Large Linewidth Source We follow the method presented in [66] and [67] for generation of a laser beam with arbitrary spectral power density. The method entails modeling the beam as radiation from a large number of independent radiators that radiate at discrete frequencies. We consider radiation from a total of N frequencies and write the n th frequency as 0n n    , where 0 is the center frequency,  is a small change 66 in frequency, and n is an integer value ranging from 2N to 2 1N  . The time dependence of the electric field for a single frequency n can be written as   12 ( ) c.c.n i t n n nE t a ib e    , (4.1) where na and nb are random variables. As a consequence of the central limit theorem, na and nb are Gaussian random variables with zero mean and variance proportional to the number of radiators at frequency n , which is determined by the power spectrum of radiation. Examples of Gaussian and Lorentzian power spectra are discussed in the following section. The time dependence of the total electric field from all frequency components is     12 ( ) c.c.n i t n n n n n E t E t a ib e      . (4.2) Setting 0n n    , the total electric field is given by      1 02 exp ( )exp c.c.n n n E t i t a ib in t      . (4.3) The term  ( )expn n n a ib in t  is a time dependent complex value. We denote this value as    t i t  :         2 1 2 exp N n n n N t i t a ib in t         . (4.4) Numerically, t is a sequence of discrete values such that t m t where m is an integer value ranging from 0 to N-1 and t is to be determined:         2 1 2 exp N n n n N m t i m t a ib inm t           . (4.5) 67 If we choose  2t N   , we get the exact definition of a discrete Fourier transform (DFT)         2 1 2 exp 2 N n n n N m t i m t a ib i nm N           . (4.6) Thus it is convenient to solve for  t and  t with the use of a DFT on the sequence  n na ib where  and  are the real imaginary part of the resulting DFT respectively. The total electric field can be written as an amplitude and phase modulated monochromatic source         1 02 exp c.c.E t t i t i t      . (4.7) Combining the time dependence described by Eq. (4.7) with the spatial geometry of the beam director, we write the electric field of a single tile as               0 0 0 1 , , ,2 0 , ,0, rect rect c.c.y x l m l m l m ik y mb ik x lb i t E x y t t i t x lb y mb E e e e                 , (4.8) where  rect  is the rectangular function defined as 0 if 2 rect( ) 1 2 if 2 1 if 2 a a a            . (4.9) The electric field of the beam director is obtained by summing over all tiles:               0 0 0 1 1 1 , ,2 1 1 0 0 , ,0, rect rect c.c.y x l m l m l m l m ik y mb ik x lb i t E x y t t i t x lb y mb E e e E e                      . (4.10) 68 Random intensity and phase fluctuations For a continuous power spectrum of radiation,  2  , the variance of na and nb are      2Var Varn n na b    . Here we consider a Gaussian power spectrum, 2 G , and a Lorentzian power spectrum, 2 L , defined as     2 02 2 1 exp 22 G GG                (4.11a)     2 2 2 0 1 L L L               (4.11b) respectively, where G and L are parameters that determine the width of each power spectrum. The power spectra are normalized such that    2 2 1G Ld d             (4.12) The linewidth,  , is defined as the full width at half maximum (FWHM) of the power spectrum. The FWHM of a Gaussian function is defined as   1 2 2 2ln 2G G   and the FWHM of a Lorentzian function is defined as 2L L   . For comparison of the two spectra, it is convenient to pick G and L such that G L       . For example, 9100 10 rad secG   and   1 2 92ln 2 117.7 10 rad secL G     which gives the relation 9235.5 10 rad secG L         . The coherence time is defined as 2ct   . Strictly speaking, this relation is specifically for a Lorentzian power spectrum. The linewidth ratio for this example is 4 0 1.24 10     . Figure 4.2 displays a comparison of two power spectra for these parameters. 69 Figure 4.2: Comparison of a Gaussian and Lorentzian power spectrum when the linewidths are matched, i.e. G L       . The Gaussian power spectrum contains 98% of the total power within  of 0 . The Lorentzian power spectrum contains 70% of the total power within  of 0 . The Lorentzian power spectrum decays more slowly than the Gaussian power spectrum and thus a greater fraction of the total power is contained at frequencies far from 0 . For example, the Gaussian power spectrum contains 98% of the total power within  of 0 , while the Lorentzian power spectrum contains 70% of the total power within  of 0 . As we will illustrate, this difference has a significant impact on the random fluctuations in the intensity and phase. We examine the behavior of the intensity and phase for a single tile at 0z  and on time scales comparable to the coherence time, ct . The intensity of a tile is defined as        0 2 2 2 2 ( ) 4 8t c c I t E t t t          , (4.13) 70 where the angled brackets indicate a temporal average performed over a time scale of 02  . It is convenient to define a normalized intensity as       2 2 0I t t t   . Due to the statistical properties of  and  we have the relation 0 ( ) 1 ct t I t . The phase of a tile is defined as      arctant t t      . (4.14) We note that strictly speaking, the arctan function has a range of  2, 2  . Numerically, we utilize the atan2 function so that the phase of a tile is defined on the interval  ,  . Figure 4.3: Intensity (a) and phase (b) fluctuations for a single tile with a Gaussian power spectrum. The intensity has random fluctuations when observed on time scales comparable to the coherence time. We begin by considering a tile with a Gaussian power spectrum defined by Eq. (4.11a). Figure 4.3a shows the intensity fluctuations at 0z  for a single tile with a Gaussian power spectrum. The intensity fluctuates randomly, but still appears as a smooth function in time. Major fluctuations occur over intervals greater than the 71 coherence time. Calculating the long term average intensity for this specific instance of  and  gives the result 0 750( ) 1.01ct t I t  . Figure 4.3b shows the phase fluctuations at 0z  for a single tile with a Gaussian power spectrum. Similar to the intensity, the phase fluctuates randomly when observed on a time scale comparable to the coherence time. Figure 4.4: Intensity (a) and phase (b) fluctuations for a single tile with a Lorentzian power spectrum. The intensity has random fluctuations when observed on time scales comparable to the coherence time. The intensity and phase as a function of time is characteristically different when the tile has a Lorentzian power spectrum described by Eq. (4.11a). Figure 4.4a shows the intensity fluctuations at 0z  for a single tile with a Lorentzian power spectrum. Similar to the case of the Gaussian power spectrum, large scale fluctuations occur over times of a few coherence time intervals. However the intensity fluctuations do not appear as a smooth function. Regardless, performing a long time average gives the result 0 750 ( ) .98 ct t I t  , which is close to the expected value of 1. 72 Figure 4.4b shows the phase fluctuations at 0z  for a single tile with a Lorentzian power spectrum. The phase displays the same random fluctuations over time scales comparable to the coherence time. Coherently combining tiles Ideally, the phase of two coherently combined beams would fluctuate synchronously and have a phase difference of zero at any instance in time. However, in practice this is generally not the case, especially for high-power lasers. Instead, coherent combining systems result in a reduction of the root mean square (RMS) phase difference between beams. To simulate a situation where the tiles are coherently combined in the transmitter plane, we develop a method to generate a set of  ,l m t and  ,l m t in which the RMS phase difference between any two tiles can be controlled. As an illustration, we consider 2 tiles whose amplitude and phase are modulated by  1 t ,  1 t ,  2 t , and  2 t . These variables are generated through the Fourier transform of the Gaussian random variables (1) (1) (2) (2), , ,n n n na b a b where we have added the superscript to distinguish between tiles. If we assume that perfect coherent combining is implemented for all frequency components within a certain frequency bandwidth, 0n   , where  is the half width of the coherent combining instrument bandwidth, then (1) (2) n na a and (1) (2) n nb b for n  . In other words, random variables associated with frequency components contained within the coherent combining bandwidth are shared between the tiles. Conversely, random variables associated with frequency components outside of the coherent combining bandwidth are uncorrelated between the tiles. The coherent combining 73 bandwidth is determined by the time scale over which the coherent combining instruments operate. We assume that the instruments have zero error and thus can precisely phase match all frequency components within the instrument bandwidth. We define 1 2, and    as the phase of tile 1, phase of tile 2, and the phase difference between the two tiles respectively. In terms of  1 t ,  1 t ,  2 t , and  2 t we have the relations      1 1 1arctant t t      , (4.15a)      2 2 2arctant t t      , (4.15b)      1 2t t t     . (4.15c) The RMS phase difference between the two tiles is defined as:   1 2 21 RMS t t N            . (4.16) We can analytically calculate RMS for the cases of incoherently combined tiles and perfectly phase matched tiles. If the two tiles are incoherently combined, then  is a random variable with a uniform distribution on the interval  0, and 1 2 2Erms     . From the properties of a uniform distribution we know the expected value,  E 2   , and variance,   2Var 12   . Using the definition of variance,     22Var E EX X X    , we calculate 2 2E 3     and 3RMS   . On the other hand, if the two tiles are perfectly phase matched then  is uniquely zero and clearly 0RMS  . These simple situations of incoherent combining and perfect phase matching correspond to 0 and  respectively. Coherently 74 combining monochromatic sources is an example of perfect phase matching. We calculate intermediary values  RMS  for both Gaussian and Lorentzian power spectra displayed in Figure 4.5. Figure 4.5: The RMS phase difference between two tiles with a Gaussian power spectrum (blue) and Lorentzian power spectrum (red) as a function of instrument bandwidth For a Gaussian power spectrum, the RMS phase difference between the two tiles is effectively zero when frequency components within 2   of 0 are matched. However, the result differs greatly when the two tiles have a Lorentzian power spectrum. Since a higher percentage of power is contained at frequencies far from 0 ,  RMS  decays much slower. From Figure 4.5, it is clear that it becomes very difficult to achieve a small RMS phase difference between two tiles with a Lorentzian power spectrum. 75 Recalling that  is the half width of the coherent combining angular frequency bandwidth, we can use Figure 4.5 to calculate the rate at which the instruments must operate to achieve any arbitrary RMS phase. For example, to achieve coherent beam combining with an RMS phase difference of 2 6 would require instruments that operate at a rate of  2 .48 2 36 GHz   and 2 2 75 GHz   for a Gaussian and Lorentzian linewidth respectively. In the following simulations, we consider the propagation of coherently combined beams with a Lorentzian power spectrum. 4.3: Propagation in turbulent atmosphere We simulate the propagation of the tile arrangement described in Section II through the atmosphere for coherently and incoherently combined tiles. For comparison, we include the case of monochromatic, phase matched tiles.As dicussed in chapter 1, we propagate the beam by numerically solving the paraxial wave equation and atmospheric turbulence is modeled as phase screens located at discrete locations along the z-axis. We use parameters that are typical for DE applications. Specifically, the center wavelength of radiation in vacuum 0 1μm  , distance to the target 5 kmL  , tile side length 2 cma  , distance between tile centers 2.2 cmb  , bandwidth 2 37.5 GHz      , and RMS phase difference between any two coherently combined tiles 2 6RMS   . We consider the lineshape to be Lorentzian. These parameters give a fractional bandwidth of 40 1.3 10    and a coherence time of 128.5 10 secct  . We consider three levels of atmospheric turbulence, 2 0nC  , 76 2 14 2 310 mnC   , and 2 13 2 310 mnC   , corresponding to vacuum, moderate turbulence and strong turbulence respectively. Using the Rytov variance, 2 2 7 6 11/6 010.5R nC L   , as a measure of turbulence strength, these conditions correspond to 2 0R  , 2 6.3R  , and 2 63R  . We use the Fried parameter,   3/5 2 2 0 0.184 nr C L    , as a measure of the transverse coherence length. Strictly speaking, 0r corresponds to the diameter of a circular area in which the RMS phase of a plane wave is 1 radian. Vacuum To observe effects exclusively associated with the linewidth of the laser sources, we begin by considering propagation in vacuum corresponding to 2 0nC  , 2 0R  , and 0r   . Figure 4.6 displays the intensity profiles at 5 kmz  for the incoherently combined beam (a), coherently combined beam (b), and monochromatic, phase matched beam (c). The intensity profiles are averaged over a time scale much greater than the coherence time, specifically 6 ns. The incoherently combined tiles diffract independently resulting in a radially symmetric, approximately Gaussian intensity profile. The energy from the incoherently combined beam is spread over a much larger area than for the coherently combined beams. The coherently combined beam produces a far field intensity pattern similar to the monochromatic beam. However, there are considerable differences. The coherently combined beam has diffracted more than the monochromatic beam which results in more beam spreading and a lower on axis intensity. 77 Figure 4.6: Intensity profiles for propagation through vacuum at 5 kmz  for the incoherently combined beam (a), coherently combined beam (b), and phase matched monochromatic beam (c). Each tile of the incoherently combined beam diffracts independently resulting in significant spreading. The coherently combined beam has an initial RMS phase difference 2 6RMS   between any two tiles. The far-field intensity profile of the coherently combined beam resembles that of the monochromatic beam. Figure 4.7 displays the intensity profile along the transverse x-axis for the incoherent beam (blue), coherent beam (green), and monochromatic, phase matched beam (red). Intensity is normalized by the on-axis intensity of the monochromatic, phase matched beam. The on axis intensity of the coherently combined beam is approximately 6.4 times greater than the on-axis intensity of the incoherently combined beam. However, the decrease in on-axis intensity of the coherently combined beam in comparison to the monochromatic beam is a fundamental limitation caused by the random temporal characteristics of the phase associated with each individual tile. A common metric used in DE applications is power in the bucket (PIB), which describes the amount of power contained in a specific area. We calculate the PIB centered at the origin as a function of bucket radius for the time averaged intensity profiles. The results are displayed in Figure 4.8 where the PIB is normalized by the total power. Again, incoherent, coherent, and monochromatic results are denoted by 78 blue, green, and red respectively. The coherently combined beam and monochromatic beam each have a distinct central lobe with a radius of approximately 7.5 cm. At a radius of 7.5 cm, the incoherent and coherently combined beams contain 37% and 78%, respectively, of the power contained by the monochromatic beam. Unlike the coherently combined and monochromatic beam, the incoherently combined beam remains radially symmetric at radii greater than the central lobe radius. Consequently, the PIB of the incoherently combined beam increases more rapidly than coherently combined beam after a radius of 7.5 cm. At a radius of 20 cm, the contained power is approximately equal for all beams. Figure 4.7: Vacuum intensity profiles along x-axis at for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). Intensity is normalized by the on-axis intensity of the monochromatic beam. The on axis intensity of the incoherently and coherently combined beams are approximately 11% and 71% of the on-axis intensity of the monochromatic beam respectively. 79 Figure 4.8: Power in the bucket (PIB) at 5 kmz  through vacuum as a function of bucket radius for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). The coherently combined beam and monochromatic beam contain significantly more power near the axis due to the distinct central lobe. The contained power within a radius of 20 cm is approximately equal for all beams. Moderate turbulence Now we consider propagation through moderate turbulence corresponding to 2 14 2 310 mnC   , 2 6.3R  , and 0 1.76 cmr  . All other parameters remain the same. An ensemble average is performed over the intensity profiles from 100 independent instances of turbulence. The intensity profile for each instance of turbulence is averaged over many coherence time intervals, 6 ns, to account for the fluctuating intensity and phase. The average intensity profiles at 5 kmz  are displayed in Figure 4.9 for the incoherently combined beam (a), coherently combined beam (b), and phase matched, monochromatic beam (c). Through comparison of the average intensity profiles we can observe a slight advantage of coherent beam combining over incoherent beam combining. These advantages include a smaller spot size and higher 80 maximum intensity. However, the effects of coherent beam combining are nearly insignificant when compared to the impact of moderate turbulence as shown in Figure 4.10. Figure 4.9: Intensity profiles for propagation through moderate turbulence, 2 14 2 310 mnC   , at 5 kmz  for the incoherently combined beam (a), a coherently combined beam (b), and phase matched monochromatic beam (c). The coherently combined beam and monochromatic beam have smaller spot size and higher on-axis intensity than the incoherently combined beam. Figure 4.10 displays the intensity profile along the x-axis for the incoherent (blue), coherent (green), and monochromatic (red) beam profiles. Again, the intensity is normalized by the on-axis intensity of the monochromatic beam in vacuum (black). The on-axis intensity of the coherently combined beam in moderate turbulence is approximately 25% of the on-axis intensity of the coherently combined beam in vacuum. Turbulence has a profound impact on the coherently combined beam in comparison to the incoherently combined beam. The on-axis intensity of the incoherently combined beam in moderate turbulence is approximately 78% of the on- axis intensity of the incoherently combined beam in vacuum. We also find that moderate turbulence significantly impacts the advantage of PIB for coherently combined beam. Figure 4.11 displays the PIB as a function of radius at 81 5 kmz  for beams propagating through moderate turbulence. The coherently combined beam and monochromatic beam no longer have a distinct central lobe as in the case of propagation through vacuum. As a result, the PIB of each beam has a similar trend. The coherently combined beam and monochromatic beam still contain more power near the axis than the incoherently combined beam, but the difference is not as significant. The advantage of coherent combining appears minor when compared to the PIB to a phase matched, monochromatic beam in vacuum. Even in moderate turbulence, it is clear that the phase distortions caused by atmospheric turbulence are dominant over the effects of coherent beam combining. Figure 4.10: Average intensity profiles in moderate turbulence, 2 14 2 310 mnC   , along x- axis at 5 kmz  for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). Intensity is normalized by the on- axis intensity of a phase-matched, monochromatic beam propagated through vacuum. For reference, the intensity profile of a phase-matched, monochromatic beam propagated through vacuum is denoted by the black curve. The on-axis intensity of the coherently combined beam is approximately twice the on-axis intensity of the incoherently combined beam. 82 Figure 4.11: Power in the bucket (PIB) at 5 kmz  through moderate turbulence, 2 14 2 310 mnC   , as a function of bucket radius for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). For reference, the PIB of the monochromatic, phase matched beam propagated through vacuum is denoted by the dashed black curve. The coherently combined beam and monochromatic beam contain more power near the axis than the incoherently combined beam, however the difference is not as significant as in vacuum. Strong turbulence Finally, we consider propagation through strong turbulence corresponding to 2 13 2 310 mnC   , 2 63R  , and 0 .44 cmr  . All other parameters remain the same. Again, an ensemble average is performed over the intensity profiles associated with 100 independent instances of turbulence and the intensity profile for each instance of turbulence is averaged over 6 ns. The average intensity profiles at 5 kmz  are displayed in Figure 4.12 for the incoherently combined beam (a), coherently combined beam (b), and phase matched, monochromatic beam (c). All three intensity profiles are nearly indistinguishable from each other. The intensity profiles of the 83 coherently combined beam and monochromatic beam have no resemblance to the vacuum far-field pattern displayed in Figure 4.6b,c. Figure 4.12: Intensity profiles for propagation through strong turbulence, 2 13 2 310 mnC   , at 5 kmz  for the incoherently combined beam (a), coherently combined beam (b), and phase matched monochromatic beam (c). The intensity profiles are nearly indistinguishable. Figure 4.13 displays the intensity profile along the x-axis for the incoherent beam (blue), coherent beam (green), and monochromatic beam (red). Again, the intensity is normalized by the on-axis intensity of the monochromatic beam in vacuum. Atmospheric turbulence is dominant over the method of beam combining; consequently, the shape of each intensity profile is nearly identical. The main noticeable difference is that the coherently combined beam and monochromatic beam have slightly higher spatial fluctuations in the time averaged intensity. 84 Figure 4.13: Intensity profiles in strong turbulence, 2 13 2 310 mnC   , along x-axis at 5 kmz  for the incoherently combined beam (blue), coherently combined beam (green), monochromatic beam (red). Intensity is normalized by the on-axis intensity of a monochromatic beam in vacuum. The intensity profiles are nearly indistinguishable. Figure 4.14 displays the PIB as a function of radius for the incoherent beam (blue), coherent beam (green), and monochromatic beam (red). The initial condition of the phase at 0z  has very little discernible impact on the PIB at 5 kmz  in strong turbulence. The advantages of coherent beam combining that were apparent in vacuum have completely disappeared in strong turbulence. 85 Figure 4.14: Power in the bucket (PIB) at 5 kmz  through strong turbulence, 2 13 2 310 mnC   , as a function of bucket radius for the incoherently combined beam (blue), coherently combined beam (green), and phase matched monochromatic beam (red). All beams deliver approximately the same power. 4.4: Discussion on beam combining In addition to the distortions caused by turbulence, there are other atmospheric effects that can modify the beam propagation. Propagation of a high-energy laser beam results in a small fraction of the laser energy being absorbed by both the molecular and aerosol constituents of air. The absorbed energy locally heats the air and leads to a decrease in the air density which modifies the ambient refractive index 0n . The change in the refractive index is 0 0( 1)TBn n    , where 0 is the ambient mass density and  is the perturbed mass density. The refractive index variation leads to a defocusing or spreading of the laser beam known as thermal blooming [14,15]. Thermal blooming can become an issue for multi-kilometer ranges and multi- hundreds of kW power levels. The effect of thermal blooming can be estimated by the 86 magnitude of the phase advance/delay caused by a change in refractive index, i.e., 2 TBn L  . For | |TBn L  , the distortion due to thermal blooming is small. Using fluid equations to describe the heating dynamics, the perturbed steady state refractive index due to thermal blooming can be estimated [15] to be TB abs TB Wn IR V   , where abs is the absorption coefficient, 4 37.5 10 cm JTB   , WV is the wind/slew velocity, R is the laser spot size, and I is the intensity. The aerosol contribution to the overall absorption coefficient can be much larger than that of the molecular water vapor. For example, in the “water window” at wavelength 1.045 μm, 5 16 10 kmabs    , while the effective aerosol contribution can be up to two orders of magnitude larger. Depending on the parameters, the phase shift due to thermal blooming can be less than the phase shift due to turbulence. For example, for 1.045 μm  , 5 cmR  , a water vapor absorption 5 16 10 kmabs    , 2250 W cmI  , at a range of 5 kmL  and a slew of 10 m sWV  , the approximate phase delay is only 2 356 . In comparison, moderate to strong turbulence can cause multiple 2 variations in phase. The result of thermal blooming is to increase the spreading of the beam, however, depending on the parameters, this effect can be much smaller than the spreading caused by moderate to strong turbulence. We have demonstrated that the linewidth of the source places a fundamental limit on the ability to coherently combine beams. The large linewidth associated with high- power lasers must be carefully considered if they are to be coherently combined. However, coherently combining beams at the transmitter plane has no benefit for DE systems in conditions of strong turbulence. We have shown that the initial RMS phase difference between tiles has negligible impact on the quality of the beam after 87 propagation through multiple kilometers of moderate or strong atmospheric turbulence. In strong turbulence, coherent beam combining at the transmitter plane merely adds to the complexity of the transmitter. Situations of strong turbulence and multi-km distances require adaptive optics to compensate for the phase distortions caused by atmospheric turbulence. Depending on the implementation of the adaptive optics, the linewidth of high-power lasers may need to be taken into consideration. For instance, coherent combining at the target can be considered an adaptive optics method because it partially compensates for distortions caused by atmospheric turbulence. In this case, the linewidth of high-power lasers places even greater limitations on coherent combining due to the transit time of light to the target and back to the receiver. If the transit time is longer the coherence time of the lasers, then it is inconceivable to coherently combine the lasers at the target without first phase matching the beams at the transmitter. There are multiple limitations that the linewidth of a laser places on the ability to coherently combine beams. Due to the dominant effect of atmospheric turbulence, it is not effective to coherently combine lasers at the transmitter plane. Incoherent combining of lasers is a much simpler approach with comparable results in moderate to strong turbulence. A more effective approach for delivering energy to a target in moderate to strong turbulence may be in using adaptive optics solutions to compensate for the turbulent distortions. 88 Chapter 5 : Bessel and Airy Beams 5.1: Introduction In this chapter, we simulate the propagation of Bessel and Airy laser beams with wavelengths of 0 1 m  through 6.4 km of weak turbulence characterized by 2 15 2/31 10nC m    , 0 1 mm , and 0 1L m . Phase screens are applied every 70 m with the first screen applied at 35z m . The transverse simulation domain is 1.2 m by 1.2 m with 2048N  cells in each direction unless otherwise stated. The transverse scale length was varied and is discussed specifically for each beam profile in the next section. The maximum propagation distance was chosen equal to the distance at which the Rytov variance, 2 2 7/6 11/6( ) 10.5R nz C z   , is unity. Longer propagation distances are considered the regime of ‘strong’ optical turbulence [1, 4]. Specifically the Rytov variance quantifies the intensity variance, 222 / 1I I I   , of a plane wave propagating through turbulent atmosphere. For weak optical turbulence, 2 0.5R  , 2 2 I R  . To validate the simulation, the linear scaling of intensity variance with Rytov variance was confirmed up to 2 ~ 0.5R for a 100 cm width Gaussian beam in a simulation domain 4 m on a side with 2048N  . The intensity variance was calculated in a circle of radius 30 cm centered at the beam centroid. The choice of 89 beam width and circle radius ensured that the intensity variance was calculated in a region were the beam propagation resembled plane wave propagation. 5.2: Gaussian beam We begin by reviewing the propagation of Gaussian beams through vacuum and turbulence. The analysis in this section forms the basis for the next sections on Bessel and Airy beams. The Rayleigh range, 2 0 /RL w  where 0w is the initial 1/ e field width, defines the length scale over which a Gaussian beam diffracts in vacuum. In particular, the on-axis intensity of a Gaussian beam drops, 2 1 0(0, ) [1 ( / ) ]RI z I z L   with a concurrent increase in the spot size, 2 1/2 0( ) [1 ( / ) ]Rw z w z L  . As discussed above, turbulence has a significant effect on beam propagation when the transverse coherence length approaches the beam diameter. Setting 0 02r w , we find an approximate length scale over which turbulent spreading and wander will cause a drop in the on-axis intensity: 2 2 5/3 1 0 0~ 0.74( )T nL C k w  . By comparing RL and TL , the relative importance of diffractive spreading and turbulence induced modifications to the beam can be ascertained. For example, if R TL L turbulence will strongly modify the beam before significant diffraction occurs. For a qualitative view of how atmospheric turbulence modifies the beam we turn to Figure 5.1. The figure displays the intensity profile of a Gaussian beam with initial spot size 0 4.5w cm after 6.4 km of propagation in (a) vacuum, (b) turbulence, and (c) averaged over 100 instances through turbulence. The color scales are normalized to the maximum in each plot. A near 5 cm wander of the beam centroid and spreading can be observed in Figure 5.1b. Figure 5.1c shows that, when ensemble averaged, the 90 beam’s width increases more rapidly than in vacuum. These modifications are expected as 6.4RL km is larger than 3.3TL km . Figure 5.1: Transverse profiles of a Gaussian beam with initial spot size 0 4.5w cm after 6.4 km of propagation in (a) vacuum, (b) turbulence with 2 15 2/31 10nC m    , and (c) averaged over 100 instances through turbulence. A more quantitative examination and scaling with 0w is presented in Figure 5.2. Results from propagation through turbulence are represented by red, from propagation through vacuum by blue, and the ratio of Fried parameter to beam diameter, 0 02r w , is displayed in green. Figure 5.2a and Figure 5.2b display 100 run ensemble averages of the normalized on-axis intensity, 0(0, ) /I z I , as a function of propagation distance for initial spot sizes of 1.8 cm and 4.5 cm respectively. The maximum of 4.5 cm was chosen such that the Rytov variance reached unity at one Rayleigh range. For 0 1.8w cm , turbulence has minimal effect on the beam’s on- axis intensity as the ratio of Fried parameter to beam diameter remains greater than unity over several Rayleigh ranges, 15TL km while 1.0RL km . In contrast, when 0 4.5w cm , 0 02 1r w  after 3.3 km well before the beam has propagated an entire 91 Rayleigh range, 6.4RL km . As a result, the on-axis intensity is ~20% lower than the vacuum value after 6.4 km of propagation. By forming the ratio 11/3 0/T RL L w  , we see that smaller beams are susceptible to intensity loss through standard vacuum diffraction, while larger beams are susceptible to spreading and wander in turbulence. In particular, noticeable turbulent modifications should occur when ~T RL L or 0 ~ 3.8w cm . This trend is illustrated in Figure 5.2c and Figure 5.2d. Figure 5.2c displays the fractional increase in root mean square (RMS) radius, 0/rmsw w , after one Rayleigh length as a function of 0w . We define the RMS radius as 2 2( ) ( , ) ( , ) rms x y I x y dxdy w I x y dxdy     (5.1) For reference the initial value is 0(0) / 1/ 2rmsw w  . In vacuum, 0/rmsw w increases by a factor of 2 regardless of initial beam width as borne out by the blue curve. Nearly the same 2 increase is observed for smaller beams, 0 ~ 2w cm , after propagation through atmospheric turbulence. But as the initial spot size increases, the effect of turbulence on the beam size becomes sizeable: the 0 ~ 4.5w cm beam expands to ~1.5 its vacuum width. Figure 5.2d shows the on-axis intensity after one Rayleigh length of propagation normalized to its initial value, 0(0, ) /RI L I , as a function of 0w . Corresponding to the increase in RMS radius, the on-axis intensity has dropped by a factor of ~2 for all 0w in vacuum and smaller 0w in turbulence. As expected, the intensity of the larger beams has dropped significantly. However, the 92 product of on-axis intensity and RMS radius squared has increased: 2 2(0, ) ( ) (0,0) (0)R RI L w L I w , suggesting that the beam has spread more from the edges than the center. Higher moments of the intensity distribution could be used to examine this effect, but we save this for future investigations. Figure 5.2: Simulation results for a Gaussian beam. (a) and (b) Comparison of the on-axis intensity in turbulence, red, and vacuum, blue, as a function of propagation distance for initial spot sizes of 1.8 cm and 4.5 cm respectively. (c) Ratio of the RMS radius at one Rayleigh length to the initial spot size as a function of initial spot size. (d) Normalized on axis intensity at one Rayleigh length as a function of initial spot size. In (a-d) the ratio of the Fried parameter to beam diameter, green, is plotted for reference. 5.3: Bessel beam Our goal now is to explore modifications to Bessel beam propagation in atmospheric turbulence. During vacuum propagation, each ring of the Bessel beam can be considered a separate beam undergoing its own diffraction. Each ring possesses near equal energy which it transports, through diffraction, outward from its initial radius 93 and inward towards the center of the beam. The inward diffraction of the rings supplies the center of the beam with energy. This process maintains the on-axis intensity until the inward diffraction of the outer ring reaches the beam center. We can estimate the diffraction length of the Bessel beam by finding the distance at which the spot size associated with outer-ring expands to the total beam size. Noting that the Bessel function zeroes approach even spacing, we write the outer ring’s spot size as ~ / ( 1)r rw R N  where rN is the number of rings and R is the aperture radius. Setting 2 1/2[1 ( / ) ]r B RR w L L  , where 2 /R rL w  , we find 2~ / ( 1)B rL R N  . The diffraction length increases with aperture area and decreases with the number of rings. A more precise length scale over which the Bessel beam diffracts is given by 20.6 / ( 1)B rL R N  , where the on-axis intensity of the Bessel beam will reach half its initial value, 0(0, ) ~ / 2BI L I , after propagating a distance BL . In the presence of turbulence, the rings lose spatial coherence and smear together, shortening the distance over which the on axis intensity remains constant. The length scale over which turbulence will cause a drop in the on-axis intensity can be approximated as: 2 2 5/3 1 0~ 0.74( )T nL C k R  . Similar to the Gaussian beam, if B TL L turbulence will modify the beam before significant diffraction occurs. For Gaussian beams we found that the cause of on-axis intensity decay depended on the beam width. Roughly speaking, small and large beams were susceptible to standard diffraction and turbulence induced modifications respectively. Based on the ratio 1 11/3/ ( 1)T B rL L N R    , we expect Bessel beams to follow the same trend for 94 fixed rN . To examine this, we fix the number of rings at 14 and vary the initial aperture radius from 6.65 cm to 23 cm. The beams are hard-apertured between adjacent intensity rings at the zeros of the Bessel function. The maximum aperture radius was chosen such that the Rytov variance was approximately unity at a distance of BL . Figure 5.3: Transverse profiles of a 15 ring Bessel beam hard-apertured at a radius of 23 cm after 6.4 km of propagation in (a) vacuum, (b) turbulence with 2 15 2/31 10nC m    , and (c) averaged over 100 instances through turbulence. Figure 5.3 displays intensity profiles of the largest aperture, 23R cm , beam after 6.4 km of propagation in (a) vacuum, (b) turbulence, and (c) averaged over 100 instances through turbulence. The color scales are normalized to the maximum in each plot. For these parameters, 6.5BL km is much larger than 0.22TL km . Figure 5.3b shows that the concentric, closed ring structure observed in Figure 5.3a has been significantly distorted by atmospheric turbulence. The rings are no longer spatially coherent with themselves or the entire beam. After 6.4 km of propagation, the transverse coherence length is 0 1.6r cm similar to the observed scale length of intensity fluctuations. 95 In Figure 5.4 the results are demarcated as before: propagation through turbulence is represented by red, propagation through vacuum by blue, and the ratio of Fried parameter to beam diameter, 0 2r R , is displayed in green. Figure 5.4a and Figure 5.4b display 100 run ensemble averages of the normalized on-axis intensity as a function of propagation distance for initial aperture radii of 6.65 cm and 23 cm respectively. The ripples in intensity result from the initial hard-aperturing. The 6.65R cm beam undergoes almost no loss of on-axis intensity due to turbulence before diffracting: 0.55BL km is smaller than 1.7TL km . With an aperture radius of 23R cm , the diffraction length far exceeds the turbulent modification length, B TL L , and the beam’s on-axis intensity drops to 90% of its vacuum value after 6.4 km of propagation. 96 Figure 5.4: Simulation results for a Bessel beam. (a) and (b) Comparison of the on-axis intensity in turbulence, red, and vacuum, blue, as a function of propagation distance for initial aperture radii of 6.65 cm and 23 cm respectively. (c) Ratio of the RMS radius at one diffraction length, BL , to the initial aperture radius as a function of initial aperture radius. (d) Normalized on axis intensity at one diffraction length as a function of initial aperture radius. In (a-d) the ratio of the Fried parameter to beam diameter, green, is plotted for reference. Comparisons of intermediate apertures radii are provided in Figure 5.4c and Figure 5.4d. Figure 5.4c displays the RMS radius normalized to the aperture radius at a distance BL as a function of aperture radius. For reference the value at the aperture is (0) / 0.57rmsw R  . Unlike the Gaussian beam, the Bessel beam’s normalized RMS radius is equal for all initial apertures after propagating a diffraction length through turbulence. The rings undergo turbulent wander and spreading, which transports energy both away from and towards the center of the beam and fills the intensity gaps between the rings. Again we treat each ring as an individual beam with spot size ~ / ( 1)r rw R N  . The total beam will undergo significant spreading in turbulence 97 after the individual rings do so: distances 2 2 5/3 1 , 0~ 0.74( )T r n rL L C k w  with 3rN  and where the subscript r refers to the individual ring. For apertures of 6.65R cm and 23R cm , , 160T rL km and , 20T rL km respectively. While the redistribution of energy within the beam does not affect its overall spreading, it does result in the decay of on-axis intensity as displayed in Figure 5.4d. As expected, the intensity decay is most severe for the larger beams, and turbulent modifications become apparent when ~T BL L or 9R cm . Thus, as with Gaussian beams, small and large beams are susceptible to standard diffraction and turbulence induced spreading and wander respectively. Fixed aperture comparison In the previous section the effect of beam size on atmospheric propagation of fixed, 14 ring Bessel beams was examined. We now fix the beam size and power, and vary the number of rings. We chose an aperture size of 30 cm , approximating the size of a beam director for directed energy applications [7]. As before, the beams are hard- apertured between adjacent intensity rings at the zeros of the Bessel function. The top plot of Figure 5.5 shows the fractional power delivered to a 30 cm aperture at distances of 1.6 km in blue, 4.0 km in red, and 6.4 km in green as a function of rings in the beam. The results in vacuum are represented by the dashed curve and turbulence by the solid curve. The delivered powers in vacuum and turbulence are nearly identical. As discussed above, over these distances, turbulence results in energy spreading within the beam without increasing the overall RMS radius beyond diffractive spreading. In particular, the distances over which turbulence modifies the 98 beam spreading for 3rN  and 14rN  are , 4.5T rL km and , 40T rL km respectively. For all three distances, the delivered power decreases as the number of rings increases. This implies, for the situation considered here, that the most effective beam for power delivery is a Bessel beam with one zero, essentially a Gaussian beam clipped at the radius of the beam director. Figure 5.5: The power delivered to a circular aperture of radius 15 cm (top) and on axis intensity (bottom) as a function of rings in a Bessel beam initially apertured at a radius of 15 cm after propagation distances of 1.6 km, blue, 4 km, red, and 6.4 km, green. The dashed and solid lines are results from propagation in vacuum and turbulence with 2 15 2/31 10nC m    respectively. The bottom plot of Figure 5.5 shows the normalized on-axis intensity as a function of rings in the beam. The distances are the same as above. With fixed power and aperture area, the initial on-axis intensity is given by 2 0 ( / 2)( 1.2)rI N R P  . Thus 99 before significant diffraction or turbulent spreading and wander, the on-axis intensity increases with number of rings as demonstrated by the blue-dashed curve. The ratio 1/ ( 1)T B rL L N   indicates that a beam’s sensitivity to on-axis intensity loss in turbulence increases with its number of rings. As illustrated by the red and green curves, the decay in on-axis intensity relative to the vacuum values increases with rN . The result of these two effects is that the ideal beam profile for applications requiring high peak intensity depends on the distance to the target. 5.4: Airy beam We now explore modifications to Airy beam propagation in atmospheric turbulence. Like the rings of the Bessel beam, each beamlet of the Airy beam can be considered a separate beam undergoing its own diffraction. During propagation through vacuum, the beamlets interfere such that the Airy beam’s center of mass follows a straight line while its intensity maximum follows a parabolic trajectory, drifting in the transverse plane. In particular, the transverse position of the intensity peak, ,r p , evolves according to 2 2 3 , , 0 ˆr r / rp i Ak z w      where ,r i is the initial position of the peak and Aw is the Airy function scaling constant, i.e. ( / ) ( / )A AAi x w Ai y w . We can estimate the Airy beam’s diffraction length, AL , by finding the distance at which an individual beamlet expands to the total beam size. To proceed, we define 2 blN as the total number of beamlets, and D as the distance from the origin to the outer edge of the th blN beamlet along a single cartesian direction. Said differently, the function ( / )AAi x w has blN zeros with the last zero at x D . Using a typical beamlet of 100 width ~ /bl blw D N , and following the same arguments used for the Bessel beam above, we find 2~ /A blL D N  . The diffraction length increases with the transverse beam area and decreases with the number of beamlets. The more precise length scale over which the Airy beam diffracts is given by 22 /A blL D N  . In particular the peak intensity of the Airy beam reaches half its initial value at a distance AL : ( , , ) ( , ,0) 2p p A p pI x y L I x y . As with Gaussian and Bessel beams, we can compute the distance at which turbulence affects the Airy beam by equating the Fried parameter with the largest transverse dimension of the beam. Approximating the Airy beam as a square, the largest dimension is the diagonal 1/22d D . Setting 0r d we find 2 2 5/3 1 0~1.2( )T nL C k D  . We expect the Airy beams to exhibit the same decay in maximum intensity as the Gaussian and Bessel beams experienced in on-axis intensity. Furthermore the ratio 1 11/3/T A blL L N D   has a similar scaling as that of the Bessel beam, suggesting that small and large Airy beams should be susceptible to standard diffraction and turbulent spreading and wander respectively. To examine the extent that turbulence modifies different size Airy beams, we fix blN at 15 and vary d from 11.3 cm to 31.1 cm. Figure 5.6 displays intensity profiles for the largest beam, 31.1d cm , after 6.4 km of propagation in (a) vacuum, (b) turbulence, and (c) averaged over 100 instances through turbulence. The color scales are normalized to the maximum in each plot. At 0z  , the maximum intensity occurs at .012x y m   . Figure 5.6a shows that the intensity peak has drifted diagonally to .12x y m  , while the tails of the beam have extended in the opposite direction, keeping the center of mass at a fixed 101 transverse position. The distortion of the intensity profile due to turbulence is demonstrated in Figure 5.6b. The intensity peak occurs near the same location as in vacuum, but the individual beamlets have undergone wander, and some are no longer discernable having smeared together. For these parameters, the Fried parameter, 0 6r cm , is still larger than a typical beamlet width, ~1.5blw cm , while 6.4AL km is much larger than 0.38TL km . Figure 5.6: Transverse profiles of a 15 zero Airy beam with 2 22 31.3d cm cm   after 6.4 km of propagation in (a) vacuum, (b) turbulence with 2 15 2/31 10nC m    , and (c) averaged over 100 instances through turbulence. Figure 5.7 uses the same color scheme as Figure 5.2 and Figure 5.4: propagation through turbulence is represented by red, propagation through vacuum by blue, and the ratio of Fried parameter to beam diagonal length, 0 /r d , is displayed in green. Figure 5.7a and Figure 5.7b display 100 run ensemble averages of the normalized on- axis intensity as a function of propagation distance for the smallest and largest beam aperture diagonals, 11.3 cm and 31.1 cm respectively. As with the Bessel Beam, the ripples in maximum intensity are due to hard-aperturing. The smallest beam undergoes standard diffraction before being affected by turbulence: 0.85AL km 102 and 2TL km . In contrast, turbulence has a pronounced effect on the on-axis intensity of the largest beam: 6.4AL km and 0.38TL  km. Figure 5.7: Simulation results for an Airy beam. (a) and (b) Comparison of the on-axis intensity in turbulence, red, and vacuum, blue, as a function of propagation distance for an aperture diagonal of 11.3 cm and 31.1 cm respectively. (c) Ratio of the RMS radius at one diffraction length, AL , to the aperture diagonal as a function of aperture diagonal. (d) Normalized on axis intensity at one diffraction as a function of aperture diagonal. In (a-d) the ratio of the Fried parameter to beam diagonal, green, is plotted for reference. Comparisons of intermediate aperture diagonal lengths are provided in Figure 5.7c and Figure 5.7d. Figure 5.7c displays the ratio of RMS radius after one diffraction length to the initial aperture diagonal, illustrating the degree to which the Airy beams spreads. For reference, the initial value is (0) / 0.44rmsw d  . As with the Bessel beam, the Airy beam’s normalized RMS radius is equal for all initial apertures after propagating a diffraction length through turbulence. The individual beamlets undergo turbulent wander and spreading filling the intensity gaps between them without 103 increasing the overall beam size. Following the same argument presented above for Bessel beams, the total beam will undergo significant spreading in turbulence after the individual beamlets do so: distances 2 2 5/3 1 , 0~ 0.74( )T bl n blL L C k w  . For diagonals of 11.3d cm and 31.1d cm , , 115T blL km and , 21T blL km respectively. Figure 5.7d shows that the maximum intensity of the Airy beams is affected by turbulence in the same way the on-axis intensity of Gaussian and Bessel beams are affected. 5.5: Discussion on quasi-non-diffracting beams We have investigated the propagation of Gaussian, Bessel, and Airy beams through atmospheric turbulence. The beam propagation was simulated using the paraxial wave equation with turbulence-induced refractive index fluctuations included through phase screens. The extent that each beam was modified by atmospheric turbulence depended on the transverse beam size. In particular, the transverse coherence length decreases with propagation distance: large aperture beams acquire transverse phase distortions before undergoing significant diffraction; small aperture beams diffract before acquiring significant transverse phase distortions. This trend held for all three beams and manifested in a drop in peak intensity. However, the nature of turbulence- induced beam spreading differed between the Bessel and Airy beams and the Gaussian beam. The simultaneous diffraction of many rings or beamlets made the Bessel and Airy beams resistant to spreading in turbulence. The propagation distance at which the Fried parameter becomes comparable to the ring or beamlet size far surpasses the distance at which it becomes comparable to the total transverse beam size. 104 We have also examined the scaling of power and intensity delivered to a target with the number of rings in a fixed aperture area Bessel beam. Despite the Bessel beams resistance to spreading, the power delivered to an area equal to the beam director decreases as the number of rings increases. As a result, the most effective Bessel beam for delivering power is a Bessel beam with zero rings, which approximates a Gaussian beam. Finally, the optimal number of rings for delivering on-axis intensity to a target depends on the distance to the target. 105 Chapter 6 : Summary and future work 6.1: Summary The work contained in this dissertation has demonstrated a new method for adaptive optics that enables the focusing of laser beams through atmospheric turbulence onto a rough surface. In addition, we have also demonstrated the limitations and shortcomings of both coherent beam combining and the propagation of quasi-non- diffracting beams. In chapter 2 we derived a novel and general relation between the field incident on a target and the field received by a point detector. We then demonstrated that the relation is valuable in the field of adaptive optics, especially when the target is uncooperative. A simulated adaptive optics system used the relation to focus a beam onto a rough surface using only the scattered field from the target to control the transmitted beam. The advantage of our method is that it does not require a beacon beam and is general to the surface properties of the target. This makes it especially suitable for applications with uncooperative targets such as directed energy or remote sensing. In chapter 3 we compared experimental and simulation results on the enhanced backscatter phenomenon. We verified previous theoretical predictions of the enhancement factor dependence on turbulence strength. Also, we developed the tilt- shift method, which replaces the standard temporal average with a spatial average. The tilt-shift method is significant because it reduces the time delay required to observe the EBS localized intensity peak. The location of EBS peak could be used to 106 accurately find the location in the detector plane which possesses the derived reciprocal properties discussed in chapter 2. The information obtained through simulation in chapter 4 demonstrates that coherent beam combining is not beneficial in situations of strong turbulence and multi-km propagation distances. It also illustrates the random intensity and phase fluctuations that are caused by the large linewidth inherent in high-power lasers. In instances of short propagation distances or weak turbulence where coherent beam combining has advantages over incoherent beam combining, these fluctuations make it very difficult to phase match the beams. We numerically determined a relation between the bandwidth of coherent instruments and the RMS phase difference between two beams. The instrument bandwidth required to coherently combine high power fiber and slab lasers is beyond the current state-of-the-art technologies. It is possible to coherently combine low power lasers; however, in general, the purpose of beam combining is to increase the total transmitted power so this has minimal applications. Chapter 5 demonstrates that quasi-non-diffracting beams are affected by turbulence in the same way as a standard Gaussian beam. Specifically, once the transverse coherence becomes comparable to the dimension of the beam, turbulent distortions become significant. We did however observe that Bessel beams did not undergo standard beam spreading due to the individual rings smearing together. Overall, the intricate structures of Bessel and Airy beams do not resist the effects of atmospheric turbulence. 107 6.2: Future work The material described in chapter 2 has the greatest potential for future work. Time dependence could be considered numerically or theoretically to possibly mitigate the occurrence of incoherent sums adding constructively as in Figure 2.4d. Also, the system could be implemented experimentally. However, there are difficulties in conducting such in an experiment. A full scale multi-kilometer experiment would require a high-power laser to enable imaging of the scattered field. Additionally, even in an experimental implementation, the limitations of adaptive optics need to be considered. These considerations include the inability to correct moderate to strong scintillation with ‘phase-only’ corrections and the limitations of optimization algorithms (as discussed in section 2.6). In an experimental setup, it would be much more difficult to work within these limits because the turbulence parameters cannot be controlled. The results obtained on coherent beam combining and quasi-non-diffracting beams imply that these are not viable methods for applications such as directed energy or free-space optics. Hence there is no clear future work in these topics for the purpose of propagating laser beams long distances, horizontally through atmospheric turbulence. If narrow linewidth, high power lasers can be developed or master oscillator – power amplifier systems can be produced with synchronous multi- aperture phase fluctuations, coherent combining may have applications in situations of weak turbulence such as propagation at higher altitudes or vertical propagation. 108 References [1] L. C. Andrews and R. L. Phillips, Laser Beam Propagation through Random Media, Bellingham, WA: SPIE, 2005. [2] R. L. Fante, "Electormagnetic beam propagation in turbulent media," Proc. IEEE, vol. 63, no. 12, pp. 1669-1692, 1975. [3] J. A. Fleck, J. R. Morris and M. D. Feit, "Time-dependent propagation of high energy laser beams through the atmosphere," Appl. Phys., vol. 10, pp. 129-160, 1976. [4] P. L. Milonni and J. H. Eberly, Laser Physics, Hoboken, NJ: Wiley, 2010. [5] R. Frehlich, "Simulation of laser propagation in a turbulent atmopshere," Appl. Opt., vol. 39, pp. 393-397, 2000. [6] S. L. Chin, A. Talebpour, J. Yang, S. Pett, V. P. Kandidov, O. G. Kosareva and M. P. Tamarov, "Filamentation of femtosecond laser pulses in turbulent air," Appl. Phys. B, vol. 74, no. 67, 2002. [7] P. A. Sprangle, A. Ting, J. Penano, R. Fischer and B. Hafizi, "Incoherent combining and atmopsheric propagation of high-power fiber lasers for directed- energy applications," IEEE J. Quant. ELectron, vol. 45, 2009. [8] C. C. Davis and I. I. Smolyaninov, "The effect of atmospheric turbulence on bit- error rate in an on-off-keyed optical wireless system," Proc. SPIE, vol. 4489, pp. 126-137, 2002. [9] H. W. Babcock, "The possibility of compensating astronomical seeing," Publ. Astron. Soc. Pac, vol. 65, p. 229, 1953. [10] V. P. Linnik, "On the possibility of reducing the influence of atmospheric seeing on the image quality of stars," Opt. Spectrosc., vol. 3, p. 401, 1957. [11] D. L. Fried, "Optical resolution through a randomly inhomogenous medium for very long and very short exposures," J. Opt. Soc. Am., vol. 56, pp. 1372-1379, 1966. [12] D. P. Greenwood, "Bandwidth specification for adaptive optics systems," J. Opt. Soc. Am., vol. 67, no. 3, pp. 390-393, 1977. [13] C. A. Primmerman, D. V. Murphy, D. A. Page, B. G. Zollars and H. T. Barclay, "Compensation of atmopsheric optical distortion using a synthetic beacon," Letters to Nature, vol. 353, pp. 141-143, 1991. [14] J. D. Ellis, Directed-energy weapons: promise and prospects, CNAS, 2015. [15] V. P. Lukin, Amospheric Adaptive Optics, Bellingham, WA: SPIE, 1995. [16] L. Puryear, J. H. Shapiro and R. R. Parenti, "Reciprocity-enhanced optical communication through atmospehric tubulence - Part II: Communication architectures and performance," J. Opt. Commun. Netw., vol. 5, no. 8, pp. 888- 900, 2013. [17] G. D. Goodno, C. P. Asman, J. Anderegg, S. Brosnan, E. C. Cheung, D. Hammons, H. Injeyan, H. Komine, W. Long and M. McClellan, "Brightness- scaling potential of actively phase-locked solid state laser arrays," IEEE J. Sel. 109 Topics Quantum ELectron., vol. 13, no. 3, pp. 460-472, 2007. [18] P. Zhou, Z. Liu, X. Wang, Y. Ma, H. Ma, X. Xu and S. Guo, "Coherent beam combining of fiber amplifiers using stochastic parallel gradient descent algorithm and its application," IEEE J. Sel. Topics Quantum ELectron., vol. 15, no. 2, pp. 240-247, 2009. [19] T. Weyrauch, M. A. Vorontsov, G. W. Carhart, L. A. Beresnev, P. A. Rostov, E. E. Polnau and J. J. Liu, "Experimental demonstration of coherent beam combining over a 7 km propagation path," Optics Letters, vol. 36, pp. 4455- 4457, 2011. [20] G. D. Goodno, S. J. McNaught, J. E. Rothenberg, T. S. Mccomb, P. A. Thielen, M. G. Wickham and M. E. Weber, "Active phase and polarization locking of a 1.4 kW fiber amplifier," Optics Letters, vol. 35, pp. 1542-1544, 2010. [21] A. E. Siegman, "How to (maybe) measure laser beam quality," OSA TOPS, vol. 17, no. 185, 1998. [22] J. Durnin, "Exact solutions for nondiffracting beams. I. The scalar theory," J. Opt. Soc. Am. A, vol. 4, no. 4, pp. 651-654, 1987. [23] P. Sprangle and B. Hafizi, "Comment on nondiffracting beams," Phys. Rev. Lett, vol. 66, no. 6, p. 837, 1991. [24] J. Arlt and K. Dholakia, "Generation of high-order Bessel beams by use of an axicon," Optics Comm., vol. 177, pp. 297-301, 2000. [25] L. Gong, Y. X. Ren, G. S. Xue, Q. C. Wang, J. H. Zhou, M. C. Zhong, Z. Q. Wang and Y. M. Li, "Generation of nondiffracting Bessel beam using digital micromirror device," Appl. Opt., vol. 52, no. 19, 2013. [26] M. A. Bandres, I. Kaminer, M. Mills, B. M. Rodrigguez-Lara, E. Greenfield, M. Segev and D. N. Chrisodoulides, "Accelerating optical beams," Optics and Photonics News, vol. 24, no. 30, 2013. [27] G. A. Siviloglou, J. Broky, A. Gogariu and D. N. Christdoulides, "Observation of Accelerating Airy Beams," Phys. Rev. Lett, vol. 99, 2007. [28] G. A. Siviloglou and D. N. Christdoulides, "Accelerating finite energy Airy beams," Opt. Lett., vol. 32, p. 979, 2007. [29] X. Ji, H. T. Eyyuboglu, G. Ji and X. Jia, "Propagation of an Airy beam through the atmosphere," Opt. Express, vol. 32, pp. 2154-2164, 2013. [30] P. Polynkin, M. Kolesik, A. Roberts, D. Faccio, P. Di Trapani and J. V. Moloney, "Generation of extended plasma channels in air using femtosecond Bessel beams," Optics Express, vol. 16, p. 15733, 2008. [31] P. Polynkin, M. Kolesik, J. V. Moloney, G. A. Siviloglou and D. N. Christodoulides, "Curved plamsa channel generation using ultraintense airy beams," Science, vol. 324, no. 5924, pp. 229-232, 2009. [32] C. Bao-Suan and P. Ji-Xiong, "Propagation of Gauss-Bessel beams in turbulent atmopshere," Chinese Phys. B, vol. 18, no. 1033, 2009. [33] I. P. Lukin, "Coherence of a Bessel beam in a turbulent atmosphere," Atmospheric and Oceanic Optics, vol. 25, pp. 328-337, 2012. [34] Y. Gu and G. Gbur, "Scintillation of Airy beam arrays in atmopsheric 110 turbulence," Opt. Letters, vol. 35, pp. 3456-3458, 2010. [35] X. Chu, "Evolution of an Airy beam in tubulence," Opt. Letters, vol. 36, pp. 2701-2703, 2011. [36] G. L. Lamb, Introductory Apllications of Partial DIfferential Equations, Hoboken, NJ: Wiley, 1995. [37] M. A. Vorontsov, G. W. Carhart and J. C. Ricklin, "Adaptive phase-distortion correction based on parallel gradient descent optimization," Optics Letters, vol. 22, no. 12, pp. 907-909, 1997. [38] V. P. Lukin and M. I. Charnotskii, "Reciprocity principle and adaptive control of optical radiation parameters," Sov. J. Quantum Electron., vol. 12, no. 5, pp. 602- 605, 1982. [39] J. H. Shapiro and A. L. Puryear, "Reciprocity-enhanced optical communication through atmopsheric turbulence - Part I: Reciprocity proofs and Far-Field power transfer optimization," J. Opt. Commun. Netw., vol. 4, no. 12, pp. 947-954, 2012. [40] J. H. Shapiro, "Reciprocity of the turbulent atmosphere," J. Opt. Soc. Am., vol. 61, no. 4, pp. 492-495, 1971. [41] J. W. Goodman, Speckle Phenomena in Optics: Theory and Applications, Englewood, CO: Roberts, 2007. [42] H. Kadono, T. Asakura and N. Takai, "Roughness and correlation-length determination of rough-surface objects using the speckle contrast," Appl. Phys. B, vol. 44, pp. 167-173, 1987. [43] V. P. Lukin, F. Y. Kanev, V. A. Sennlkov, N. A. Makenova, V. A. Tartakovskli and P. A. Konyaev, "Phase and amplitude-phase control of a laser beam propagating in the atmopshere," Quantum Electron., vol. 34, no. 9, pp. 825-832, 2004. [44] C. A. Primmerman, T. R. Price, R. A. Humphreys, B. G. Zollars, H. T. Barclay and J. Herrmann, "Atmospheric-compensation experiments in strong-scintillation conditions," Appl. Opt., vol. 34, pp. 2081-2088, 1995. [45] M. C. Roggemann and D. J. Lee, "Two-deformable-mirror concept for correcting scintillation effects in laser beam projection thorugh the turbulent atmopshere," Appl. Opt., vol. 37, no. 21, pp. 4577-4585, 1998. [46] V. A. Banakh and V. L. Mironov, LIDAR in a Turbulent Atmosphere, Dedham, MA: Artech House, 1987. [47] G. Welch and R. Phillips, "Simulation of enhanced backscatter by a phase screen," J. Opt. Soc. Am. A, vol. 7, no. 4, pp. 578-584, 1990. [48] E. Jakeman, "Enhanced backscattering through a deep random phase screen," J. Opt. Soc. Am. A., vol. 5, no. 10, pp. 1638-1648, 1988. [49] J. H. Churnside and J. J. Wilson, "Enhanced backscatter of a reflected beam in atmospheric turbulence," Appl. Opt., vol. 32, pp. 2651-2655, 1993. [50] Y. N. Barabankov, Y. A. Kravtsov, V. D. Ozrin and A. I. Saichev, "Enhanced backscattering in optics," Progess in Optics, vol. 29, 1991. [51] R. A. Murphy and R. L. Phillips, "Reciprocal path-scattering effects for ground based, monostatic laser radar tracking a space target through turbulence," Appl. 111 Opt., vol. 36, pp. 5997-6004, 1997. [52] A. G. Vinogradov, Y. A. Kravtsov and V. I. Tatarskii, "Enhanced backscattering from bodies immersed in a random inhomogenous medium," Izv. VUZ Radiofiz., vol. 16, pp. 1064-1070, 1973. [53] A. S. Gurvich and S. S. Kashkarov, "Enahnced backscattering in a turbulent medium," Izv. VUZ Radiofiz., vol. 20, no. 5, pp. 794-796, 1977. [54] S. S. Kashkarov, "Enhancement of the average backscattered intensity in a turbulent atmosphere," Izv, VUZ Radiofiz., vol. 26, no. 1, pp. 44-48, 1983. [55] M. S. Belen'kii and V. L. Mironov, "Diffraction of optical radiation on a mirror disc in a turbulent atmosphere," Quantum Electron., vol. 5, pp. 38-45, 1972. [56] J. E. Harvey, S. P. Reddy and R. L. Phillips, "The precision pointing and tracking through random media by explointing the enhanced backscatter phenomenon," Appl. Opt., vol. 35, pp. 4220-4228, 1996. [57] G. Welch and W. T. Rhodes, "Imaging through a random screen by enhanced backscatter: a physical interpretation," Appl. Opt., vol. 32, pp. 2692-2696, 1993. [58] A. Consortini, Y. Y. Sun, C. Innocenti and Z. P. Li, "Measuring inner scale of atmospheric turbulence by angle of arrival and scintillation," Opt. Commun., vol. 216, pp. 19-23, 2003. [59] V. I. Tatarskii, Wave Propagation in a Turbulent Atmosphere, Moscow: Nauka Press, 1967. [60] R. G. Lane, A. Glindemann and J. C. Dainty, "Simulation of a Kolmogorov phase screen," Waves in Random Media, vol. 2, no. 3, pp. 209-224, 1992. [61] A. Ishimaru, "Wave propagation and scattering in random media and rough surfaces," Proceedings of the IEEE, vol. 79, no. 10, pp. 1359-1366, 1991. [62] D. H. Nelson, D. L. Walters, E. P. MacKerrow, M. J. Schmitt, C. R. Quick, W. M. Porch and R. R. Petrin, "Wave optics simulation of atmospheric turbulence and reflective speckle effects in CO2 lidar," Appl. Opt., vol. 39, no. 12, pp. 1857- 1871, 2000. [63] L. C. Andrews, R. L. Phillips, C. Y. Hopen and M. A. Al-Habash, "Theory of optical scintillation," J. Opt. Soc. Am. A, vol. 16, no. 6, p. 1417, 1999. [64] A. Kellerer, Assesing time scales of atmopsheric turbulence at observatory sites, Paris: Astro-Physics, Universite Paris-Diderot, 2007. [65] L. Mandel and E. Wolf, Optical Coherence and Quantum Optics, Cambridge: Cambridge University Press, 1996. [66] N. Ruggieri, D. Cummings and G. Lachs, "Simulation of superposed coherent and chaotic radiation of arbitrarty spectral shape," J. Appl. Phys., vol. 43, no. 3, 1972. [67] G. Vannucci and M. C. Teich, "Computer simulation of superposed coherent and chaotic radiation," Appl. Optics, vol. 19, no. 4, pp. 548-553, 1980. [68] E. R. Mendez, M. A. Ponce, V. Ruiz-Cortes and Zu-Han, "COherent effects in the scattering of light from random surfaces with symmetry," Opt. Lett., vol. 16, no. 3, pp. 123-125, 1991. 112 [69] M. I. Charnotskii, "Common omissions and misconceptions of wave propagation in turbulence: discussion," J. Opt. Soc.Am A, vol. 29, no. 5, pp. 711-721, 2012. [70] J. D. Jackson, Classical Electrodynamics, 3rd Edition, Hoboken, NJ: Wiley, 1999. [71] C. C. Davis, Lasers and Electro-Optics, 2nd Edition, New York: Cambridge University Press, 2014. [72] J. W. Goodman, Introduction to Fourier Optics, Englewood, CO: Roberts & Co., 2005.