Abstract:
An array of acoustic signal sources and sensors are positioned in a fluid. Acoustic signals propagated through the medium are subjected to phase differences as they pass through fluid currents. An eikonal description of the acoustic waves propagating through the fluid medium is used to derive the fluid currents. Fluid currents are measured in a multiplicity of acoustic paths extending from the acoustic sources to the sensors, using acoustic phase tomography.

Description:
BACKGROUND OF THE INVENTION 
     1. Field of the Invention 
     This invention relates to a method and apparatus for measuring fluid currents, particularly water currents, in order to obtain surface or bulk current distributions. More particularly, this invention relates to the use of accoustic transmissions in fluids, with eikonal approximations and tomographic techniques used to derive current measurement from the phase information of acoustic signals. 
     2. Description of the Related Art 
     There exists a need for on-site determination of water currents. Knowledge of water currents is important in monitoring and predicting conditions in various marine ecosystems, including fisheries. Such knowledge can also be beneficial in determining the transport of sediment and pollutants. Knowledge of pollutant transport can help to mitigate damage to fisheries and other wildlife habitats. 
     Previous methods for determining currents have included flow meters and high-frequency radar. Flow meters can provide on-site measurements, but not at the water&#39;s surface where wind effects are dominant. High-frequency radar can be used to measure surface currents, but the range resolution of such radar does not provide details at scales less than a kilometer. 
     Acoustic waves have been used to detect and map underwater objects such as submarines and ocean floor profiles. It is also known that a sound wave traveling in a fluid contains information about physical properties of the fluid itself. Previous acoustic tomographic techniques have used an inversion method or inverse calculation from the travel time of the acoustic wave to determine properties of the fluid such as temperature, salinity, and density. Other techniques have concentrated on the intensity and amplitude of the acoustic wave as sources of information. 
     SUMMARY OF THE INVENTION 
     The invention determines surface or bulk currents using phase information contained in acoustic wave propagated through the fluid. An array of acoustic sources and acoustic sensors is used to obtain a tomographic image of currents in the fluid. 
     When a sound wave passes between two points through a fluid which has no movement relative to the transmission and receiving points, the sound wave is received at the receiving point with a certain phase for a particular time. However, if there are fluid currents along the acoustic trajectory between the two points, the phase of the acoustic signal can be shifted. Using an eikonal description of acoustic waves propagated through the fluid, relevant information about the current distributions in the fluid may be derived from the phase information of the received acoustic waves. 
     In a preferred embodiment, an array of coherent underwater acoustic transmitters and receivers are used. The transmitters will typically operate in the 250 to 3000 Hz frequency band, although other frequencies may also be used. Where multiple transmitters will be transmitting simultaneously, the transmitters may each operate at a different frequency. This allows signals from various transmitters to be readily distinguished from one another. 
     When surface currents are the main object of interest, these transmitter and receiver arrays can be located just beneath the water&#39;s surface in an essentially 2-dimensional array. Where bulk currents are of interest, the transmitters and receivers can be located at various depths to form a 3-dimensional array. 
     The invention can be used in various liquids, but has particular application in bodies of water such as oceans, rivers, lakes and estuaries. 
    
    
     DESCRIPTION OF THE DRAWINGS 
     FIG. 1a is an elevation view of a simple transmitter and receiver apparatus useful for explaining the invention&#39;s derivation of current information from acoustic phase information; 
     FIG. 1b is an elevation view of a simple transmitter and receiver apparatus, similar to that shown in FIG. 1a; 
     FIG. 2 is a top plan view of one embodiment of the invention comprising a 2-dimensional array of transmitters and receivers; 
     FIG. 3 is a perspective view of another embodiment of the invention comprising a 2-dimensional array of transmitters and receivers positioned beneath the surface of a body of water; and 
     FIG. 4 is a perspective view of a further embodiment of the invention which includes a 3-dimensional array of transmitters and receivers for determining bulk currents. 
    
    
     DETAILED DESCRIPTION OF THE INVENTION 
     FIGS. 1a and lb show a transmitter T and a receiver R positioned in a fluid beneath the fluid surface 10. An acoustic signal 12 is propagated through the fluid from the transmitter T to the receiver R. In FIG. 1a there is no fluid current along the unperturbed trajectory 14, so that the acoustic signal 12 is received at an expected phase angle φ 0 , which in the example of FIG. 1a is π/2 or 90° for the time T 0 . This expected (or zero-current) phase angle φ 0  is a function of the distance s between the transmitter T and receiver R, the phase and frequency of the acoustic signal as transmitted, and the velocity of sound for the fluid. 
     However, where there is a fluid current along the unperturbed trajectory, the received acoustic signal is received with a phase angle which is shifted from the so-called zero-current phase angle. FIG. 1b shows the transmitter T and receiver R positioned beneath the fluid surface 10, with the fluid having a current 16 along the unperturbed trajectory 14. The passage of the acoustic signal through the areas of current along the unperturbed trajectory is impeded or assisted by the current, depending upon whether the current is in the direction of or against the acoustic signal. In FIG. 1b the current is against the acoustic signal, i.e., the current is in the direction away from the receiver and toward the transmitter. The net velocity of the acoustic signal is thus retarded. 
     When the acoustic signal 12 is received at the receiver it has a phase angle φ 1 , in this case π or 180°, which is shifted from the zero-current phase angle φ 0 . For the acoustic signal of FIG. 1b, the phase angle is shifted by +90° from the phase angle for the zero-current example of FIG. 1a. Analysis of this phase information gives the average projection of the velocity along the acoustic trajectory. 
     For a single transmitter and receiver, as shown in FIG. 1a and lb, it can be difficult to determine the true amount of the phase shift. Because the phase angle can only be measured through one cycle or oscillation, i.e., through 360°, there is the chance that a phase shift may actually be greater than one cycle, or that the phase shift might be negative. For example, a phase shift measured as 90° might actually be a phase shift of 450° or -270°. 
     To avoid any confusion about the actual versus measured phase shift, the acoustic signal could be a pulsed signal, and any difference in the pulsed signal&#39;s travel time (as compared to the zero-current travel time) could be used to determine the actual phase shift. Another alternative is to have the transmitter transmit at several different frequencies and to use the cross-correlated measured phase shifts of all frequencies to determine the actual phase shifts. Another approach makes use of multiple transmitters and receivers with multiple and intersecting acoustic paths; by cross-correlating the information along each acoustic path, the actual phase shifts for each acoustic signal are determined. 
     Once the actual phase shift is known, the average velocity along the acoustic trajectory between a single transmitter T and single receiver R can be determined. By increasing the number of transmitters and receivers, the current distribution (as opposed to merely the average velocity) along each of several acoustic trajectories can be determined. Through the use of known tomographic techniques, the current distribution within a 2-dimensional or 3-dimensional area can be determined from the phase shifts. 
     Assuming that the currents change over distances which are large relative to the wavelengths of the acoustic signals that are used to probe them, the effect of the currents on the phases of the acoustic waves can be determined by means of the adiabatic, or eikonal, approximation. In this approximation: ##EQU1## where ω is the radian frequency of the acoustic wave and is related to k by a known dispersion relation 
     
         ω=ω(k,r)                                       (4) 
    
      denotes a vector, ∇ is the gradient function, ti is time, i is the complex number symbol, k is the wave vector of the acoustic wave, dr/dt denotes the group velocity of the acoustic signal, n(r,t) is the density of the acoustic medium, and k(r) and k=|k(r) | is determined by solving the dispersion relation. 
     For acoustic waves in a fluid, the dispersion relation is simply 
     
         ω=kV.sub.s                                           (5) 
    
     where V s  is the speed of sound. If (a) the fluid is in motion with a position-dependent velocity v(r,t), (b) the distance over which v(r,t) varies is large compared to the acoustic signal length, and (c) the time over which v(r,t) varies is long compared to ω -1 , then the dispersion relation of equation 5 can be applied in the frame of reference moving with the fluid as follows: 
     
         ω-k·V(r,t)=kV.sub.s.                        (6) 
    
     Then we have 
     
         ∇ω=∇(k·V(r,t))=(k·∇)V(r,t)+kx(∇xV(r,t))                                   (7) 
    
     where x is the vector product symbol and ∇x is the curl operation, and ##EQU2## To the first order in the velocity field v(r,t), this says that the phase occurring in equation (3), evaluated along an acoustical signal path, is simply ##EQU3## where k 0  =ω/V s , s is the distance along the unperturbed trajectory, and V k .sbsb.0 (f,t) is the component of v(r,t) in the direction of the unperturbed trajectory. Accordingly, in the eikonal approximation, the necessary information about the current flow is contained in the phase of the acoustic signal received by the detector when a transmitter transmits an acoustic wave through the region of interest. Using an eikonal approximation for phase information reconstruction, the current flow along the acoustic trajectory can be determined. 
     With several acoustic sources placed at different locations, the current flow distribution in space can be determined by standard tomographic techniques. Tomographic algorithms permit the extraction of position-dependent quantities from integrals of these quantities along trajectories. Roughly, if there are M sources and N detectors, then MxN position-dependent quantities can be determined in this manner. So for properly placed arrays of acoustic arrays and detectors, the three Cartesian components of the current vector can be determined by (M×N)/3 positions. 
     Equation (9) shows that the phase change along an acoustical signal path is approximately ##EQU4## Thus, if for a 2-dimensional case the region is discretized into an N×M grid, and N×M paths are chosen through the region, there will be N×M values of the phase change along the different paths. The NxM equations giving the phase change in terms of the V k .sbsb.0 at the N×M values of  defining the grid coordinates can then be solved for the N×M values of V k .sbsb.0 at the grid coordinates by conventional tomographic techniques, which yield information about the spatial dependence of a quantity in terms of integrals of the quantity. One example of a suitable tomographic algorithm, which employs a matrix inversion for this purpose, is described in Somerstein et al., &#34;Radio frequency geotomography for remotely probing the interiors of operating mini- and commercial-sized oil-shake retorts&#34; Geophysics, Vol. 49, No. 8, Aug. 1984, pages 1288-1300. 
     The use of computer aided tomography (CAT) in medical, geophysical and other applications is well known. One such technique, useful when phase information is of interest, is embodied in this simple equation, which relates the received power at a designated receiver, P R , to the transmitted power T  at a designated transmitter: ##EQU5## where Γ 2  =efficiency of the transmitter/receiver pair 
     γ=α+iβ= propagation vector 
     α= attenuation coefficient 
     β=phase shift 
     Υ=distance between transmitter and receiver K T  (⊖ T ) =(amplitude) gain pattern for transmitter in the direction ⊖ T  away from its beam axis. 
     K R  (⊖ R )=(amplitude) gain pattern for receiver in the direction ⊖ R  away from its beam axis. 
     Equation (10) is equally applicable to the description of acoustic waves as well as electromagnetic waves. 
     Discretization of equation (10) yields: ##EQU6## The index i, enumerates the grid points encountered by the discretized path, and the index k enumerates different paths between the transmitter and receiver. These paths are treated as zero-width beams between the transmitting and receiving devices. Equation (11) is equivalent to a system of linear algebraic equations: 
     
         A x=b, 
    
     and the input dataset b is available to measurement. If there were a unique path from transmitter to receiver, the structure matrix A could be obtained and inverted by convergent iterative methods. In the present case, however, this is not true, because hydrophonic devices do not emit narrow beams. The received signal should therefore be regarded as resulting from propagation along several paths between receiver and transmitter. 
     For clarity, we continue to restrict the description of the process to two spatial coordinates. Then the angular coordinate k is one-dimensional, and the signal detected at the ith receiver from the jth transmitter is obtained as a convolution over the index k. 
     Equation (11) in this case is of the form: ##EQU7## where 
     
         P.sub.(i,j) =P.sub.Rc δ.sup.2 K.sub.T.sub.2 (⊖.sub.R)K.sub.R.sup.2 (⊖.sub.R)          (14) 
    
     Here g(n) is a window function introduced to minimize spurious effects due to discretization. One example of such a window function is the Hamming window defined by: 
     
         g(m)=0.54-0.46 cos {2πn/(k-1)} for m=0,1 . . . k-1      (15) 
    
     The limits of summation in equation (14) are extended to plus or minus infinity, because the window function vanishes outside the beam width. ##EQU8## where P(i,j,k)=p(i,j)exp {-ξ∂(i,j,(e)(17)α(i,j} Because equation (13) is of convolution type, the solution can be effected very efficiently by the Discrete Fourier Transform method. The result is: 
     
         P(z)=P.sub.RC (Z)/.sub.G (Z)                               (18) ##EQU9## where and N is the dataset size P(Z), P.sub.RC (Z) and G (Z) are FFT&#39;s of, P(on index k), P.sub.RC (i,j) and g(k) respectively. 
    
     By reversing the Fourier Transform we can now obtain a space sequence: 
     
         P(i,j,k)=P(i,j)exp {-Σ&#39;Υ(ijk)αΥ(i,j)}(19) 
    
     which gives the signal through a single path. Index i labels the receiver j labels the transmitter, and k labels the path. 
     To complete the solution, replace P with p in equation (11), which then describes a zero width path. In this case the &#34;algebraic iterative reconstructive technique&#34; such as described by the following articles; R. Gordon., 1974, entitled &#34;A tutorial on ART (algebraic reconstruction techniques): Trans. Inst. of Elect. &amp; Electron. Engineers, Nucl., Sci., NS-21, 78-93, G. T. Herman, A. Lent and S. W. Rowland, 1973 entitled, ART: Mathematics and Applications (A report on the mathematical foundations and of the applicability to real data of the algebraic reconstruction techniques), Journal of Theoretical Biology, Vol. 42, pages 1 -32, and the article to D. L. Lager and R. J. Lytle, entitled &#34;Determining a subsurface electromagnetic profile from high frequency measurements by applying reconstruction technique algorithms: Radio Science, Vol. 12, pages 249-260, is directly applicable and produces best-fit estimates of the propagation vector Υ as a function of spatial coordinates, along that particular path. Since there are multiple receivers and multiple transmitters, the M×N paths can be chosen in such a way as to span the whole region with any reasonable resolution cell. 
     The foregoing discussion of the proposed method is applicable to any kind of wave propagation in an inhomogeneous medium. To relate the result to propagation of acoustic waves, simply refer back to equation (9) where it is clear that the total phase shift: ##EQU10## is to be identified with the imaginary part of γdr in equation (10). Attenuation of acoustic waves (i.e., the real part of γ) can also be obtained as a by-product. 
     In applications where the sources operate simultaneously, each source may be operated at a different frequency from the other sources. In this way, signals from different sources may be distinguished from each other, so that there is no ambiguity about the origin of a signal received at a detector. Preferred frequencies are in the band extending from about 250 to 3000 Hz. However, other frequencies may also be used, depending upon the particular application. 
     A precise knowledge of the distance s between the transmitters T and receivers R is important. As one option, all transmitters and receivers can be firmly anchored in position, as where the transmitters and receivers are secured to fixed platforms such as dock pilings or rigidlysecured buoys. Another option is to employ a system to precisely track the positions of all transmitters and receivers, and adjust the values of s to compensate for any changes in position. 
     FIG. 2 shows one embodiment of the invention. The apparatus 20 shown comprises a 2-dimensional array of receivers R and transmitters T. The array includes a line 22 of receivers R and a line 24 of transmitters T. The acoustic paths 26 run from each transmitter to each receiver. Thus, the apparatus 20 covers an area loosely defined at its opposite edges by the transmitters T and receivers R. The outer regions of the area near the transmitters and receivers have a lower density of acoustic paths 26 than the central region. To ensure that sufficient information is gathered, the array preferably extends beyond the boundaries of the region of significant interest. 
     A 2-dimensional array such as that shown in FIG. 2 will determine velocity distributions only in the plane of the array. That is, velocity distributions will be determined in only two dimensions. Given M sources and N detectors, M×N position-dependent quantities can be determined. For 2-dimensional arrays, the two planar components of the current vector can be determined at (M×N)/2 positions. Thus, for the 6×6 array shown in FIG. 2, the two planar components of the current vector can be determined at 18 positions. 
     FIG. 3 shows another embodiment, similar to that shown in FIG. 2, in which the system 30 comprises an essentially 2-dimensional array of transmitters T and receivers R located just below the water&#39;s surface 32. Such an array is useful for determining surface currents. The receivers R and transmitters T are shown suspended underneath respective buoys 34 and 35 by respective lines 36 and 37. The receivers R may be of a variety of types, so long as they can operate in an underwater environment. Hydrophones are effective receivers for the underwater applications of this invention. 
     The determination of the phase differentials between the transmitted and received signals can be performed at any convenient location, such as on one of the transmitter or receiver buoys, on a separate buoy or floating platform, or on land if close enough. There is a similar latitude in locating a site for calculating the fluid currents once the phase differentials are known, although typically both functions will be performed at the same location. A conventional IQ comparator can be used to determine the phase differentials, and a standard computer for the current calculation. In FIG. 3, signals corresponding to the transmitted signals are delivered from the transmitters T to their respective buoys 35, such as by a fiber optic linkage, while the received signals are delivered from the receivers R to their buoys 34. The optical fibers from the transmitted T and receivers R are bundled in respective fiber optic cables 38 and 39, which deliver the signals to an IQ comparator 40 or other mechanism for determining the phase differentials. The phase differential output from IQ comparator 40 is delivered to a computer microprocessor 41 that is programmed to perform the tomographic algorithm that calculates the fluid currents. This information can then be transmitted by another fiber optic cable (not shown) to a desired pickup location. Communications channels other than fiber optics, such as RF broadcasts from the transmitter and receiver buoys to a computation station and from the computation station to a separate monitoring location, could also be employed. 
     The array as shown includes four receivers R and four transmitters T. However, this is merely one example of a potential array; the specific numbers and positions of the transmitters and receivers is dependent upon the particular application. 
     FIG. 4 shows a system 42 which comprises a 3-dimensional array such as may be used in determining bulk currents. In the embodiment shown, a series of transmitters T are suspended on a line 43 at varying depths beneath each floating transmission buoy 44. Similarly, a series of receivers R are suspended on a line 46 beneath each floating receiver buoy 48. A communications link between both the transmitters and the receivers and a computation station, such as a fiber optic network similar to that illustrated in FIG. 3, would also be included. 
     The acoustic tranjectories 50 between transmitters and receivers cover a 3-dimensional region that is roughly enclosed at the water&#39;s surface 32 by the buoys 44, 48. In the illustration of FIG. 4, in which nine total transmitters and nine total receivers are used, the three Cartesian components of the current vector can be determined at (9×9)/3, or 27, positions. 
     It should be noted that essentially any number and combination of transmitters and receivers may be used. The choice of numbers and positions for the transmitters and receivers depends upon the specific application. For example an array used for monitoring currents in an estuary will normally have a different physical layout than an array used for monitoring near-shore ocean currents. In general, with increased numbers of transmitters and receivers, more precise current distribution information may be obtained. 
     The above described embodiments are intended to illustrate the principles of the invention, but not to limit its scope. Other embodiments and variations will be apparent to those skilled in the art, and may be made without departing from the spirit and scope of the invention as defined in the following claims.