Exisitng synthetic aperture radar image formation systems and techniques are complex and computationally demanding. Range migration algorithm (RMA) methods require a straight and level flight path of the aircraft. Although aircraft path compensation to the first order is commonly built into the algorithm, image distortion may still be present at both ends of the range swath when the flight path deviates from certain bound. A second issue associated with RMA is that motion compensation becomes very computational demanding for a SAR system employing a wide beam in both elevation and azimuth, such as in ground penetration radar.
Polar format (PF) algorithm methods can tolerate more flight path turbulence by performing motion compensation repeatedly on a large number of image patches. However, since the wave front for the radar pulses follows an arc instead of a straight line its tomography-type processing leads to geometric distortion and phase distortion. Consequently, pixel phase continuity only presents within each small image patch. This greatly reduces its capability in coherent change detection (CCD). In addition, a close range SAR requires many image patches in Polar Format which makes the processing very inefficient.
In radar signal processing, both RMA and PF perform pulse-to-pulse processing during the initial data collection phase. The processing along the azimuth dimension is held until pulses equivalent to the full synthetic aperture have been acquired. This tends to leave a much higher computation load during batch processing. Consequently, it requires a processor with a relatively high computation throughput rate. For RMA, it requires a large memory buffer to store the range compressed data. Then, there is a need to perform a data corner turning process in order to access data in the azimuth direction. These demand a large size memory and fast I/O bandwidth from the processor.
The technology described herein relates to synthetic aperture radar (SAR) image formation systems and techniques that produce high quality images for all major performance matrices. One embodiment is a synthetic aperture radar (SAR) imaging method. The method includes applying range pulse compression to SAR in-phase and quadrature (IQ) data for each image block of a SAR image grid to generate range-compressed IQ data for each image block. The method also includes applying pulse range frequency (PRF) decimation to the range-compressed IQ data to generate PRF-decimated range-compressed IQ data for each image block. The method also includes applying motion compensation to the PRF-decimated range-compressed IQ data to generate motion-compensated data for each image block. The method also includes computing first stage image values at image grid point intersections of iso-range lines and vertical grid lines for each image bock based on the motion-compensated data for each image block. The method also includes computing second stage image values for the image grid point intersections by interpolation using the first stage image values at the image grid point intersections. The method also includes correcting image phase of the second stage image values for the image grid point intersections to generate phase-corrected image values for each image block. The method also includes generating a full-resolution SAR image by summing the phase-corrected image values for each image block.
In some embodiments, applying motion compensation includes performing a fast Fourier transform in the SAR image grid range dimension. In some embodiments, the fast Fourier transform in the SAR image grid range dimension is performed in accordance with:
where S(k,i) is a two-dimensional range compressed SAR response with k being range bin index and i being pulse index, Np is number of pulses in a processing pulse group, and Nfft is the least integer power of 2 greater than the number of slant range bins covering the processed image block.
In some embodiments, the method includes calculating:
S
1(n,i)=S′(n,i)·Ref(n)
by multiplying S′(n,i) with a reference function of the form:
Ref(n)==exp(−i2π·n1·ΔR(i)/Δrs),
where n=[1, 2, 3 . . . , Nfft], n1=[0, 1, 2 . . . , Nfft/2−1, −Nfft/2, −Nfft/2+1, . . . , −1], ΔR(i)=R(i)−R(ic), where R(i) is distance between a reference scatter and the SAR at pulse i which is within a range of [1, 2, 3, . . . , Np], and ic=Np/2+1.
In some embodiments, the motion-compensated data for each image block is determined by performing a fast Fourier transform in the range dimension in accordance with:
In some embodiments, the first stage image values are computed in accordance with:
where i′=[−Np/2, . . . , −1, 0, 2, . . . Np/2−1], fd(Rj,n) is Doppler frequency at (Rj,n), Rj=∥{right arrow over (x)}ac(ic)−({right arrow over (x)}b+n·{circumflex over (p)}az+β·{circumflex over (p)}rg)∥, {right arrow over (x)}b is position of each image block center, {right arrow over (x)}ac(i) is position of the SAR aircraft at the i-th pulse in the pulse group, {circumflex over (p)}az is azimuth unit vector of the image grid, {circumflex over (p)}rg is range unit vector of the image grid, n is azimuth coordinate index, ic=Np/2+1, Np is number of pulses in a processing pulse group, where n is the along-track pixel index and β is cross-track pixel index associated with the pixel having a slant range of Rj.
In some embodiments, the second stage image values are computed in accordance with:
where R(m,n) is distance between grid point (m,n) and {right arrow over (x)}ac(ic), and wk(Rj,m) are the interpolation coefficients determined based on the range of the point (m,n) and Rj.
In some embodiments, the full-resolution SAR image is generated in accordance with:
where Ri
Another embodiment is a synthetic aperture radar (SAR) imaging system. The system includes one or more processors, and a memory. The memory includes executable code representing instructions that when executed cause the system to apply range pulse compression to SAR in-phase and quadrature (IQ) data for each image block of a SAR image grid to generate range-compressed IQ data for each image block. The memory includes executable code representing instructions that when executed cause the system to apply pulse range frequency (PRF) decimation to the range-compressed IQ data to generate PRF-decimated range-compressed IQ data for each image block. The memory includes executable code representing instructions that when executed cause the system to apply motion compensation to the PRF-decimated range-compressed IQ data to generate motion-compensated data for each image block. The memory includes executable code representing instructions that when executed cause the system to compute first stage image values at image grid point intersections of iso-range lines and vertical grid lines for each image bock based on the motion-compensated data for each image block. The memory includes executable code representing instructions that when executed cause the system to compute second stage image values for the image grid point intersections by interpolation using the first stage image values at the image grid point intersections. The memory includes executable code representing instructions that when executed cause the system to correct image phase of the second stage image values for the image grid point intersections to generate phase-corrected image values for each image block. The memory includes executable code representing instructions that when executed cause the system to generate a full-resolution SAR image by summing the phase-corrected image values for each image block.
In some embodiments, the memory includes executable code representing instructions that when executed cause the system to apply motion compensation by performing a fast Fourier transform in the SAR image grid range dimension.
In some embodiments, the fast Fourier transform in the SAR image grid range dimension is performed in accordance with:
where S(k,i) is a two-dimensional range compressed SAR response with k being range bin index and i being pulse index, Np is number of pulses in a processing pulse group, and Nfft is the least integer power of 2 greater than the number of slant range bins covering the processed image block.
In some embodiments, the method includes calculating:
S
1(n,i)=S′(n,i)·Ref(n)
by multiplying S′(n,i) with a reference function of the form:
Ref(n)==exp(−i2π·n1·ΔR(i)/Δrs),
where n=[1, 2, 3 . . . , Nfft], n1=[0, 1, 2 . . . , Nfft/2−1, −Nfft/2, −Nfft/2+1, . . . , −1], ΔR(i)=R(i)−R(ic), where R(i) is distance between a reference scatter and the SAR at pulse i which is within a range of [1, 2, 3, . . . , Np], and ic=Np/2+1.
In some embodiments, the motion-compensated data for each image block is determined by performing a fast Fourier transform in the range dimension in accordance with:
In some embodiments, the first stage image values are computed in accordance with:
where i′=[−Np/2, . . . , −1, 0, 2, . . . Np/2−1], fd(Rj,n) is Doppler frequency at (Rj,n) Rj=∥{right arrow over (x)}ac(ic)−({right arrow over (x)}b+n·{circumflex over (p)}az+β·{circumflex over (p)}rg)∥, {right arrow over (x)}b is position of each image block center, {right arrow over (x)}ac(i) is position of the SAR aircraft at the i-th pulse in the pulse group, {circumflex over (p)}az is azimuth unit vector of the image grid, {circumflex over (p)}rg is range unit vector of the image grid, n is azimuth coordinate index, ic=Np/2+1, Np is number of pulses in a processing pulse group, where n is the along-track pixel index and β is cross-track pixel index associated with the pixel having a slant range of Rj.
In some embodiments, the second stage image values are computed in accordance with:
where R(m,n) is distance between grid point (m,n) and {right arrow over (x)}ac(ic), and wk(Rj,m) are the interpolation coefficients determined based on the range of the point (m,n) and Rj.
In some embodiments, the full-resolution SAR image is generated in accordance with:
where Ri
The SAR methods and systems described herein (hereinafter “technology”) can provide one or more of the following advantages. One advantage of the technology is that the technology preserves the phase of all ground scatters within an image and synthesizes high resolution images from lower resolution images from the image domain. The technology provides high image quality in terms of impulse response and geometric accuracy. The technology is insensitive to flight path turbulence. The technology provides phase continuity within the entire image. The technology is reduces the computation load of the batch processing by performing both azimuth processing and range processing during pulse collection phase. The technology is breaks the SAR processing task into coarse grain jobs such that it may fully utilize the processing power of a processor that has a very large number of processing nodes.
The foregoing features of various embodiments of the invention will be more readily understood by reference to the following detailed descriptions in the accompanying drawings.
An antenna of the aircraft 112 illuminates a patch 114 as the aircraft moves along the flight path 108. The antenna illuminates the patch 114 by transmitting radar pulses towards the patch 114 and by receiving corresponding radar return pulses to generate the SAR image array. Typically a complete array has 2n pulses (e.g., 256, 512, 2048 pulses) so it is compatible with Fast Fourier Transform (FFT) processing to be performed on the pulses to generate the SAR image.
Table 1 includes variables used in describing the embodiment. In this image processing technique, the position of a frame reference point at the top-left corner 134 of the image grid 132 is determined in accordance with:
{right arrow over (x)}
FM
ref
={right arrow over (x)}
CRP
−L
FM/2·{right arrow over (p)}rg−LFM/2·{right arrow over (p)}az EQN. 1,
and {right arrow over (x)}m,n is the position of an image grid 132 point with grid index of m and n, which is determined in accordance with:
{right arrow over (x)}
m,n
={right arrow over (x)}
FM
ref+(m·{right arrow over (p)}rg+n·{right arrow over (p)}az)·Δx EQN. 2.
The Doppler of {right arrow over (x)}m,n with respect to {right arrow over (x)}ac(i) is rdi,m,n, and is given by:
The controller 240 includes a communication module 216, one or more input devices 220, one or more output devices 228, one or more display devices 224, one or more processors 232, and memory 236. The modules and devices described herein can, for example, utilize the processor 232 to execute computer executable instructions and/or the modules and devices described herein can, for example, include their own processor to execute computer executable instructions. It should be understood the controller 240 can include, for example, other modules, devices, and/or processors known in the art and/or varieties of the described modules, devices, and/or processors.
The communication module 216 includes circuitry and code corresponding to computer instructions that enable the computing device to send/receive signals to/from the antenna 204. For example, the communication module 216 provides commands from the processor 232 to the transmitter 208 to control how the antenna 204 transmits radar pulses during operation. The communication module 216 also, for example, receives data corresponding to the radar return pulses received by the receiver 212. The received data can be, for example, stored by the memory 236 or otherwise processed by the processor 232.
The input devices 220 receive information from a user (not shown) and/or another computing system (not shown). The input devices 220 can include, for example, a keyboard, a scanner, a microphone, a stylus, a touch sensitive pad or display. The output devices 228 output information associated with the controller 240 (e.g., information to a printer, information to a speaker, information to a display, for example, graphical representations of information). The processor 232 executes the operating system and/or any other computer executable instructions for the controller 240 (e.g., executes applications). The memory 236 stores a variety of information/data, including profiles used by the controller 240 to specify how the system 200 generates the full SAR image. The memory 236 can include, for example, long-term storage, such as a hard drive, a tape storage device, or flash memory; short-term storage, such as a random access memory, or a graphics memory; and/or any other type of computer readable storage.
The method then includes applying range pulse compression (step 308) to the IQ data for each image block of the SAR image grid to generate range-compressed IQ data for each image block. In collecting the SAR IQ data, a long-duration frequency moduled pulse is typically transmitted. The pulse has the property that, when filtered with a matched filter, the result is a narrow pulse in which all the pulse energy has been collected to the peak value. Matched filtering of the received pulses is called range compression. Range compression is performed on each range line of SAR IQ data using a Fast Fourier Transform (FFT) method. The frequency domain range matched filter is generated once, and it is applied to each range line. The range matched filter is typically computed or generated from a replica of the transmitted pulse. In addition, the range matched filter frequency response typically includes an amplitude weighting to control sidelobes in the range impulse response.
The method then includes applying pulse range frequency (PRF) decimation (step 312) to reduce the pulse sampling rate. PRF decimation is accomplished first by convolving the range-compressed IQ data with the kernel of an along-track band-pass filter to remove the spectra outside of the imaged area. It is then followed by a down sampling step to reduce the pulse rate.
The method then includes applying motion compensation (step 316) to the PRF-decimated range-compressed IQ data to generate motion-compensated data for each image block. The range-compressed radar return IQ data (from a point-like scatter) has an impulse-like function with its peak position located at a delay determined by its distance to the aircraft position. In this embodiment, to simplify the description, range distance is used instead of delay time and the radar range gate begins when the pulse is leaving the aircraft. At the peak of the impulse response, there is also a phase term given by the two-way distance and the center wavelength. With respect to the pulse index, this phase profile mainly includes a constant term which is the phase at the center pulse, a linear term proportional to the Doppler of the scatter, and a quadratic term.
The motion compensation process use a reference scatter in the image block center to determine the range and phase profile to be compensated. It performs range migration correction such that scatter peak responses are aligned at the same range bin as that of the center pulse. This range migration correction applies to all scatters in the image block. In addition, the quadratic phase is also removed. Implementation of motion compensation also takes into account the slight variation of the quadratic phase as a function of range. The range migration correction is implemented by FFT of the range compressed data block, associated with the image block, reference multiply, and inverse FFT of the composite spectra.
In this implementation, S(k,i) represents the two-dimensional range compressed SAR response with k being the range bin index and i being the pulse index. The reference function generated for motion compensation is given by:
Ref(n)==exp(−i2π·n1·ΔR(i)/Δrs), EQN. 4
where n=[1, 2, 3 . . . , Nfft], n1=[0, 1, 2 . . . , Nfft/2−1, −Nfft/2, −Nfft/2+1, . . . , −1], ΔR(i)=R(i)−R(ic), where R(i) is the distance between the reference scatter and the aircraft at pulse i which is within the range of [1, 2, 3, . . . , Np], ic is the index of center pulse (i.e., ic=Np/2+1), N is the number of pulses in a processing pulse group, and Nfft is the least integer power of 2 greater than the number of slant range bins covering the processed image block.
In motion compensation, the first step is to perform a fast Fourier transform for S(k,i) in the range dimension, in accordance with:
Then, this signal spectrum is multiplied with the reference function in accordance with:
S
1(n,i)=S′(n,i)·Ref(n) EQN. 6
The last step is to perform an inverse fast Fourier transform in the range dimension in accordance with:
Within the processed image block, the range migration of scatters at different azimuth location will not be completely removed since their range migration profiles are slightly different from that of the center scatter. However, if the residual range migration is less than a small fraction of a range bin, it has very little impact to the final impulse response. However, if the chosen image block size is large, this residual range migration will be significant such that it will degrade the final 2D impulse response. In such case, the solution is to introduce the keystone process [x] after motion compensation to completely remove the range migration of all scatters in the image block. Keystone processing removes the linear term of target range migration in the SAR data. This is performed in a group of SAR data pulses after range filtering by rescaling the slow time following in accordance with:
t1=fc/(f+fc)·t, EQN. 8
where t is the slow time, fc is the center frequency and f+fc is the frequency of the rescaled data; and as further described in “SAR Imaging of Moving Targets”, Perry, et al. IEEE Transactions on Aerospace and Electronic Systems, Vol. 35, No. 1, January 1999; the entire contents of which is hereby incorporated by reference.
The method then includes computing first stage image values (step 320) at image grid point intersections of iso-range lines and vertical grid lines for each image bock based on the motion-compensated data for each image block.
Computing the first stage image values involves compute image values at the pixel locations 408. The ‘range’ in the iso-range term is the distance of the pixel location 408 to the antenna (e.g., antenna the radar 116 of
As shown in
R
j
=∥{right arrow over (x)}
ac(ic)−({right arrow over (x)}b+n·{circumflex over (p)}az+β·{circumflex over (p)}rg)∥, EQN. 9
where {right arrow over (x)}b is the location of the block center, n is the along-track pixel index and β is cross-track pixel index associated with the pixel having a slant range of Rj. After the location of the circular dot is obtained, we compute its Doppler frequency using EQN. 3. Let the azimuth coordinate index be n, the Doppler frequency at (Rj,n) is expressed as fd(Rj,n).
The image value of this circular point is obtained by filtering the motion compensated data at the given Doppler frequency for the range bin corresponding to Rj, where Rstart≦Rj≦Rend This filtering process preserves the phase, at the center pulse, of the scatter. The 1D motion compensated data at the range of Rj is determined in accordance with Smocomp(i), i=[1, 2, . . . , Np]. The filtering process is expressed in accordance with:
where i′=[−Np/2, . . . , −1, 0, 2, . . . Np/2−1]. The use of instead of i is for preserving the scatter phase. The notation of the output Sstage
The next stage of the method includes computing second stage image values (step 324) for the image grid point intersections by interpolation using the first stage image values at the image grid point intersections. The is accomplished by interpolation along the first dimension of Sstage
where R(m,n) is distance between grid point (m,n) and {right arrow over (x)}ac(ic), and wk(Rj,m) are the interpolation coefficients determined based on the range of the point (m,n) and Rj.
Before storing the computed image value into the image buffer, a phase correction step is performed. The method includes correcting image phase of the second stage image values (step 328) for the image grid point intersections to generate phase-corrected image values for each image block. This is performed because the same scatter processed by different pulse groups has different phase values as shown in
In this step of process, the phases of all image grid points within the image block is computed by their distance to the aircraft position at the center of its own pulse group. This process creates ramp like phases, as shown in
R
i
(l),m,n
=∥{right arrow over (x)}
ac(ic(l))−{right arrow over (x)}m,n∥ EQN. 12
The full resolution image is formed by coherently summing 332 these phase corrected responses to generate a sharply defined scatter image in accordance with:
producing a high resolution response as illustrated as the output of the summing block in
In some embodiments, the method includes the optional steps of autofocus (step 336), impulse response shaping (step 340), detection/gain compensation (step 344), and outputting the SAR image data (step 348). The autofocus process (step 336) estimate the nonlinear phase associated with targets in the image and then removes the phase error to yield a well focused image. Various autofocus processes can be used, including those described in EP 1367409B1, issued to Cho, “System and Method for Focusing SAR images,” U.S. Published Application No. 2003/0222807 by Cho, “Focusing SAR Images Formed by RMA with Arbitrary Orientation,” U.S. Pat. No. 6,781,541, issued to Cho, “Estimation and Correction of Phase for Focusing Search Mode SAR Images Formed by Range Migration Algorithm,” U.S. Pat. No. 8,009,079, issued to Connell et al., “Methods For Two-Dimensional Autofocus in High Resolution Radar Systems,” and U.S. Published Application No. 2006/0109164 by Cho et al., “Efficient Autofocus Method for Swath SAR;” The entire contents of each of these references hereby incorporated by reference.
The impulse response shaping process (step 340) can be performed by, for example, applying a Taylor weighted taper in both range and azimuth dimension. Impulse response shaping is often applied to achieve a more uniform resolution over the entire image area by applying weighting coefficients to the pulses. In each dimension, this process is performed in accordance with
S(x)=FFF−1(FFT(S(x))·WTaylor(X)) EQN. 14
The detection process of step 344 can be performed by, for example, converting the complex SAR image pixel in-phase and quadrature data (I,Q) into real SAR image pixel value S in accordance with:
S=√{square root over (I2+Q2)} EQN. 15
The gain compensation process of step 344 can be performed by, for example, removing the gain modulation due to the antenna pattern onto the SAR image based on the knowledge of the antenna pointing direction and the antenna pattern. This process can be performed in accordance with:
S
1(x,y)=S(x,y)·G−1(x,y) EQN. 16
The above-described systems and methods can be implemented in digital electronic circuitry, in computer hardware, firmware, and/or software. The implementation can be as a computer program product that is tangibly embodied in non-transitory memory device. The implementation can, for example, be in a machine-readable storage device and/or in a propagated signal, for execution by, or to control the operation of, data processing apparatus. The implementation can, for example, be a programmable processor, a computer, and/or multiple computers.
A computer program can be written in any form of programming language, including compiled and/or interpreted languages, and the computer program can be deployed in any form, including as a stand-alone program or as a subroutine, element, and/or other unit suitable for use in a computing environment. A computer program can be deployed to be executed on one computer or on multiple computers at one site.
Method steps can be performed by one or more programmable processors, or one or more servers that include one or more processors, that execute a computer program to perform functions of the disclosure by operating on input data and generating output. Method steps can also be performed by, and an apparatus can be implemented as, special purpose logic circuitry. The circuitry can, for example, be a FPGA (field programmable gate array) and/or an ASIC (application-specific integrated circuit). Modules, subroutines, and software agents can refer to portions of the computer program, the processor, the special circuitry, software, and/or hardware that implement that functionality.
Processors suitable for the execution of a computer program include, by way of example, both general and special purpose microprocessors, and any one or more processors of any kind of digital computer. Generally, a processor receives instructions and data from a read-only memory or a random access memory or both. The essential elements of a computer are a processor for executing instructions and one or more memory devices for storing instructions and data. Generally, a computer can be operatively coupled to receive data from and/or transfer data to one or more mass storage devices for storing data. Magnetic disks, magneto-optical disks, or optical disks are examples of such storage devices.
Data transmission and instructions can occur over a communications network. Information carriers suitable for embodying computer program instructions and data include all forms of non-volatile memory, including by way of example semiconductor memory devices. The information carriers can, for example, be EPROM, EEPROM, flash memory devices, magnetic disks, internal hard disks, removable disks, magneto-optical disks, CD-ROM, and/or DVD-ROM disks. The processor and the memory can be supplemented by, and/or incorporated in special purpose logic circuitry.
Comprise, include, and/or plural forms of each are open ended and include the listed parts and can include additional parts that are not listed. And/or is open ended and includes one or more of the listed parts and combinations of the listed parts.
One skilled in the art will realize the invention may be embodied in other specific forms without departing from the spirit or essential characteristics thereof. The foregoing embodiments are therefore to be considered in all respects illustrative rather than limiting of the invention described herein. Scope of the invention is thus indicated by the appended claims, rather than by the foregoing description, and all changes that come within the meaning and range of equivalency of the claims are therefore intended to be embraced therein.