WO2022070242A1 - 画像解析装置および画像解析方法 - Google Patents

画像解析装置および画像解析方法 Download PDF

Info

Publication number
WO2022070242A1
WO2022070242A1 PCT/JP2020/036843 JP2020036843W WO2022070242A1 WO 2022070242 A1 WO2022070242 A1 WO 2022070242A1 JP 2020036843 W JP2020036843 W JP 2020036843W WO 2022070242 A1 WO2022070242 A1 WO 2022070242A1
Authority
WO
WIPO (PCT)
Prior art keywords
value
image analysis
low
unit
pixel
Prior art date
Application number
PCT/JP2020/036843
Other languages
English (en)
French (fr)
Inventor
大地 田中
Original Assignee
日本電気株式会社
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 日本電気株式会社 filed Critical 日本電気株式会社
Priority to JP2022553249A priority Critical patent/JP7416277B2/ja
Priority to PCT/JP2020/036843 priority patent/WO2022070242A1/ja
Priority to US18/027,687 priority patent/US20230377201A1/en
Publication of WO2022070242A1 publication Critical patent/WO2022070242A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/97—Determining parameters from multiple pictures
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G01S13/9023—SAR image post-processing techniques combined with interferometric techniques
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration

Definitions

  • the present invention relates to an image analysis device and an image analysis method applicable to an image analysis system that performs image analysis based on the phase of the received electromagnetic wave of a synthetic aperture radar.
  • Synthetic Aperture Radar (SAR) technology is a technology for transmitting and receiving electromagnetic waves while moving objects such as artificial satellites and aircraft, and obtaining observation images equivalent to antennas with large openings.
  • Synthetic aperture radar is used, for example, to process reflected waves from the ground surface and analyze altitude, ground surface displacement, and the like.
  • the analysis device inputs a time-series SAR image (SAR data) obtained by a synthetic aperture radar, and analyzes the input SAR image in a time-series manner.
  • Interference SAR analysis is an effective method for analyzing altitude and surface displacement.
  • the phase difference of the radio wave signals constituting a plurality of (for example, two) SAR images taken at different times is calculated. Then, based on the phase difference, the change in the distance between the flying object and the ground that occurs during the shooting period is detected.
  • Patent Document 1 describes an analysis method using a coherence matrix.
  • the coherence matrix represents the correlation of pixels at the same position in multiple complex images.
  • Coherence is calculated by the complex correlation of pixels at the same position in a plurality of SAR images out of S (S ⁇ 2) SAR images.
  • the pair of SAR images be (p, q) and the components of the coherence matrix be cp , q .
  • p and q are values of S or less, respectively, and indicate any of S SAR images.
  • the phases ⁇ p and q are calculated.
  • the absolute value of the value obtained by averaging exp ( ⁇ j ⁇ p, q ) for a plurality of pixels in a predetermined region including the pixel to be calculated for coherence is the component c p, q of the coherence matrix.
  • the coherence matrix contains information such as dispersion that can estimate the degree of phase noise.
  • the phases ⁇ p and q correlate with the displacement speed and the difference in shooting time.
  • the displacement is estimated based on the mean value of the phase difference. It is possible to verify the accuracy of the displacement analysis by using the amount of phase noise. Therefore, the coherence matrix can be used for displacement analysis.
  • the phases ⁇ p, q correlate with the elevation of the object to be analyzed and the distance between the projectiles (eg, the distance between two imaging positions of the projectile).
  • the altitude is estimated based on the average value of the phase difference. It is possible to verify the accuracy of elevation analysis by using the amount of phase noise. Therefore, the coherence matrix can be used for elevation analysis.
  • Patent Document 1 describes a method of fitting a model such as displacement to a coherence matrix and restoring a phase excluding the influence of noise. Further, Non-Patent Document 1 also discloses a similar method.
  • noise included in the phase difference can be reduced.
  • the noise included in the phase difference is as small as possible, so that more noise is desired to be removed.
  • An object of the present invention is to provide an image analysis device and an image analysis method capable of increasing the degree of phase noise reduction.
  • the image analysis apparatus is a pixel selection means for selecting pixels at a plurality of positions in one image among a plurality of images recorded in the same area, and a complex as an evaluation value when the evaluation function is optimized.
  • a low-dimensional compression means that low-dimensionally compresses a vector into a low-dimensional space
  • an expansion means that returns the compression result by the low-dimensional compression means to the original pixel space and calculates a spatial correlation phase estimate, and an evaluation value of the spatial correlation phase. It includes an optimization means that optimizes the evaluation function by bringing the estimated value closer to the pixel at the position selected by the pixel selection means.
  • pixels at a plurality of positions in one image among a plurality of images recorded in the same area are selected, and a complex vector as an evaluation value when the evaluation function is optimized is low-dimensional.
  • Low-dimensional compression into space the low-dimensional compressed result is returned to the original pixel space, the spatial correlation phase estimate is calculated, and the evaluation value is brought closer to the spatial correlation phase estimate and the pixel at the selected position. By doing so, the evaluation function is optimized.
  • the image analysis program is a process of selecting pixels at a plurality of positions in one image of a plurality of images recorded in the same area by a computer, and as an evaluation value when an evaluation function is optimized.
  • the process of optimizing the evaluation function is executed by bringing the image closer to the pixel at the selected position.
  • the degree of phase noise reduction can be further increased.
  • the image analyzer has a first evaluation function (observation signal evaluation function) based on a weighted observation signal (observation pixel value) and a noise-free pixel value, and a phase difference in the image.
  • a second evaluation function spatial correlation evaluation function
  • the image analysis device calculates the phase estimated to be noise-free by using an evaluation function that integrates the observation signal evaluation function and the spatial correlation evaluation function.
  • “Integration” is, for example, "sum”.
  • evaluation function it means an evaluation function in which an observation signal evaluation function and a spatial correlation evaluation function are integrated.
  • the image analysis device of this embodiment uses an inverse matrix of the absolute values of the coherence matrix with a negative sign as the weight.
  • the inverse matrix of the absolute value of the coherence matrix with a negative sign is an example, and the image analyzer may use other parameters as weights.
  • the image analysis device may include a coherence matrix calculation unit.
  • the coherence matrix calculation unit calculates the coherence matrix as follows, for example.
  • the coherence matrix calculation unit calculates the coherence matrix C for, for example, S sheets of SAR images (complex image: including amplitude and phase information) stored in the SAR image storage unit. For example, let the pair of SAR images be (p, q) and the components of the coherence matrix C be c p, q . p and q are values of S or less, respectively, and indicate any of S SAR images.
  • the coherence matrix calculation unit calculates the phases ⁇ p and q (specifically, the phase difference) for the pair of SAR images.
  • the coherence matrix calculation unit uses the values obtained by averaging exp ( ⁇ j ⁇ p, q ) of a plurality of pixels in a predetermined region including the pixel to be calculated for coherence as the components c p, q of the coherence matrix C.
  • a p , A q , and exp ( ⁇ j ⁇ p, q ) may be averaged, where the intensity in the SAR image p is Ap and the intensity in the SAR image q is A q .
  • each element of the matrix may be divided by the value obtained by dividing the sum of the diagonal components by S. ..
  • the matrix obtained as the average of A p , A q , and exp (-j ⁇ p, q ) is multiplied by a diagonal matrix having the diagonal component -1 / 2 as a diagonal component from the left and right. You may.
  • the coherence matrix calculation unit can calculate the coherence matrix by various methods.
  • the observation signal evaluation function any of the functions that can be expressed in various ways can be adopted.
  • the observation signal evaluation function for example, has a weight W (in the following embodiments, a coherence matrix) in addition to at least a noise-free pixel value x, that is, an estimated pixel value x, as an evaluation value.
  • W in the following embodiments, a coherence matrix
  • the image analysis device may include a spatial correlation prediction unit.
  • the spatial correlation prediction unit predicts the correlation coefficient of pixels in an image and generates a spatial correlation matrix K which is a matrix having the correlation coefficient as an element.
  • the spatial correlation matrix K is a matrix having a pixel correlation coefficient as an element.
  • the spatial correlation matrix K is an example of the representation of the spatial correlation, and other representations may be used.
  • the spatial correlation prediction unit creates (calculates) spatial correlation (correlation of pixels in the analysis target region) based on prior information (known data) regarding the analysis target region (may be the entire region) in the SAR image. ..
  • Prior information is given to the spatial correlation prediction unit in advance. As the prior information, for example, the following information is exemplified.
  • -Rough shape of the analysis target for example, the structure to be the displacement analysis
  • the analysis target for example, the structure to be the displacement analysis
  • the image analysis system to which the image analysis device 10 is applied.
  • a spatial correlation is created so that the pixels have a smooth correlation in the stretching direction of the structure.
  • -Weights in graphs nodes: pixels, edges (weights): correlations
  • a spatial correlation matrix K having a value based on the weight as an element is created.
  • the spatial correlation evaluation function including the spatial correlation matrix K generated based on the prior information corresponds to a function indicating that the phase difference is smoothly distributed in the image.
  • prior information is not limited to the above example, and may be other types of information as long as the spatial correlation of the analysis target can be predicted.
  • the spatial correlation evaluation function any of the functions that can be expressed in various ways can be adopted.
  • the spatial correlation evaluation function is, for example, at least the spatial correlation matrix K corresponding to a kind of weight (actually, the inverse matrix of the spatial correlation matrix K) in addition to the pixel value x as the evaluation value. It is an evaluation function (for example, a product-sum operation formula) having K -1 ) as a parameter.
  • the image analysis device is selected from all the pixels (pixels 1 to N) in the predetermined region, instead of executing the optimization process using the observation signal evaluation valence.
  • the optimization process using the observation signal evaluation valence is executed for the M (M ⁇ N) pixels.
  • the evaluation function includes an observation signal evaluation function and a spatial correlation evaluation function
  • the image analysis device decomposes the spatial correlation matrix K in the spatial correlation evaluation function into low dimensions (low rank approximation). Execute the optimization process. More specifically, the image analysis apparatus substantially decomposes the inverse matrix K -1 of the spatial correlation matrix K into low dimensions.
  • the image analyzer can perform low-dimensional decomposition using, for example, a technique called Nystrom approximation. In that case, the image analyzer decomposes the inverse matrix K -1 of the spatial correlation matrix K in a low dimension using the equation (1).
  • G is a matrix of d ⁇ d (d ⁇ N).
  • ⁇ is an N ⁇ N diagonal matrix.
  • FIG. 1 is a block diagram showing a configuration example of the image analysis device of the first embodiment.
  • the image analysis device 10 shown in FIG. 1 includes a pixel selection unit 110, a low-dimensional compression unit 120, a decompression unit 130, and an optimization unit 140.
  • the image analysis device 10 may include a coherence matrix calculation unit and a spatial correlation prediction unit (not shown in FIG. 1).
  • the SAR image is stored in the SAR image storage unit 200 in advance.
  • the weights used in the optimization process are stored in advance in the weight matrix storage unit 300.
  • the weight matrix as weights is stored in the weight matrix storage unit 300.
  • the SAR image storage unit 200 and the weight matrix storage unit 300 may be installed outside the image analysis device 10 or may be included in the image analysis device 10.
  • the pixel selection unit 110 selects a plurality of pixels in a certain SAR image from a plurality of SAR images in which the same area is recorded. Then, the pixel selection unit 110 identifies a pixel at the same position as the pixel in the above-mentioned one SAR image in another SAR image. The pixel selection unit 110 designates the position of the specified pixel.
  • the low-dimensional compression unit 120 maps a pixel at a position designated by the pixel selection unit 110 to a low-dimensional space.
  • the expansion unit 130 maps the pixel at the position designated by the pixel selection unit 110 from the low-dimensional space to the original pixel space.
  • the optimization unit 140 uses a weight matrix for a pair of the above-mentioned one SAR image and another SAR image.
  • the pixel selection unit 110 designates a common position in one SAR image and another SAR image. Then, the optimization unit 140 has a phase so as to be close to both the calculation result output from the expansion unit 130 and the pixel value (observation pixel value) of the observation pixel at the position designated by the pixel selection unit 110. Pixels (represented by complex vectors; corresponding to x samp, s described later) are derived.
  • the optimization unit 140 optimizes (for example, maximizes) the spatial correlation evaluation function to bring the complex vector having the phase closer to the calculation result output from the expansion unit 130. Further, the optimization unit 140 optimizes (for example, maximizes) the observation signal evaluation function to bring the complex vector closer to the pixel value of the observation pixel at the position designated by the pixel selection unit 110. Also, when the observation signal evaluation function is optimized, the complex vector approaches the pixel value of the observation pixel in the region where the noise is relatively small.
  • the pixel selection unit 110 specifies the position of the pixel, but the corresponding pixel (specifically, the pixel value) is read from the SAR image storage unit 200, and the read pixel is referred to as a low dimension. It may be supplied to the compression unit 120, the expansion unit 130, and the optimization unit 140.
  • the weight W is stored in advance in the weight matrix storage unit 300. Further, the inverse matrix K -1 of the spatial correlation matrix K represented by the equation (1) is set in advance in the image analysis device 10.
  • the low-dimensional compression unit 120 uses the matrix GV T which is a part of the inverse matrix K -1 of the spatial correlation matrix K.
  • the expansion unit 130 uses the matrix V which is a part of the inverse matrix K -1 of the spatial correlation matrix K. Therefore, the matrix GV T and the matrix V may be set in the image analysis device 10.
  • the pixel selection unit 110 generates M random numbers (M ⁇ N) (step S101).
  • N indicates the total number of pixels of the SAR image.
  • M is a value of 2 or more.
  • the pixel number of the selected pixel is nm ( m : 1 to M).
  • nm becomes any one of 1 to M.
  • S indicates the total number of SAR images. Let the image number be s (s: 1 to S). Further, nm may be expressed as samp .
  • the pixel selection unit 110 uses each of the generated random numbers as the pixel number of the selected pixel. That is, the pixel selection unit 110 randomly selects pixels.
  • the expansion unit 130 multiplies the previous f s by V samp (step S102).
  • the previous f s means the operation result of the process of step S104 executed immediately before.
  • Each of the V samps ( samp : 1 to M) is a row vector generated by extracting n 1 row, n 2 rows, ..., N M rows from the matrix V. That is, the expansion unit 130 calculates V samp f s for each pixel of the row corresponding to the pixel numbers n 1 , n 2 , ..., N M by using the corresponding row vector in the matrix V.
  • V samp f s may be referred to as a spatial correlation phase estimate.
  • V samp f s corresponds to the value (complex number data) obtained by multiplying VGV T in Eq. (1) by x samp, s . Since "s" in "samp, s" is an image number, for example, "samp, 1" is an index that specifies the selected pixel in the SAR image of image number 1.
  • the optimization unit 140 optimizes (for example, maximizes) the evaluation function (step S103).
  • the evaluation function includes an observation signal evaluation function including a weight W and an observation pixel value y as parameters in addition to a pixel value x that is assumed to have no noise as an evaluation value, and a spatial correlation matrix K.
  • it is an evaluation function integrated with the spatial correlation evaluation function having the inverse matrix K -1 of the spatial correlation matrix K) as a parameter.
  • the optimization unit 140 converts x samp, s into W samp (in the present embodiment, the inverse matrix of S ⁇ S coherence matrices corresponding to the selected pixels) with respect to the observation signal evaluation function.
  • W samp in the present embodiment, the inverse matrix of S ⁇ S coherence matrices corresponding to the selected pixels
  • the low-dimensional compression unit 120 compresses the x samp, s optimized by the optimization unit 140 using the GV T samp (step S104). That is, the low-dimensional compression unit 120 multiplies x samp, s and GV T samp .
  • the low-dimensional compression unit 120 further multiplies the multiplication result by (N / M).
  • the value obtained by multiplying the multiplication result by (N / M) is expressed as f s . Since the processes of steps S102 to S104 are executed for M pixels selected from all N pixels, the multiplication result of x samp, s and GV T samp is the original value (M). / N) It is doubled.
  • the low-dimensional compression unit 120 multiplies the multiplication result by (N / M) and returns it to the original value.
  • the low-dimensional compression unit 120 has a built-in memory.
  • the low-dimensional compression unit 120 is configured to temporarily store fs in a memory (not shown in FIG. 1) and pass it to the expansion unit 130.
  • step S106 If the end condition is satisfied, the process ends (step S106). If the end condition is not satisfied, the process returns to step S101.
  • the pixel selection unit 110 selects a pixel group (which may be partially overlapped) different from the previously selected pixel group in the process of step S101 executed next.
  • the end condition is, for example, that V samp f s as the spatial correlation phase estimated value output from the expansion unit 130 is less than a predetermined value predetermined with respect to the value obtained in the previous process of step S102. That is. That is, the end condition is that the fluctuation amount of V samp f s is small. In other words, the termination condition is that it is determined that V samp f s has converged to the optimum value. Further, the end condition may be that x samp, s output from the optimization unit 140 is less than a predetermined value predetermined with respect to the value obtained in the previous process of step S103. That is, the end condition may be that the fluctuation amount of x samp, s is small.
  • the end condition may be determined that x samp, s has converged to the optimum value. It should be noted that such an end condition is an example, and another condition, for example, a condition that the processes of S101 to S104 are executed a predetermined time may be used.
  • the low-dimensional compression unit 120 calculates GV T samp x samp, s .
  • the expansion unit 130 expands a value ( fs ) obtained by multiplying GV T samp x samp, s by (N / M) by V.
  • the optimization unit 140 optimizes the spatial correlation evaluation function using a pixel group different from the previously selected pixel group with respect to the observation signal evaluation function. By repeating the processing by the low-dimensional compression unit 120, the decompression unit 130, and the optimization unit 140, x samp, s converges to the optimum value.
  • the image analysis device 10 obtains a phase close to the phase of the observed signal by optimizing the spatial correlation evaluation function to smoothly distribute the phase difference and optimizing the observed signal evaluation function. Therefore, the degree of phase noise reduction can be further increased.
  • the amount of memory used can be reduced as compared with the case where the space correlation matrix K is not decomposed in a low dimension.
  • the amount of memory used is (the number of dimensions of f s ⁇ the number of images) + the number of selected pixels ⁇ the number of images 2 ).
  • the optimization unit 140 executes the optimization process for the pixels randomly selected from all the pixels with respect to the observation signal evaluation function. Therefore, the amount of calculation of the optimization unit 140 is reduced.
  • the pixel selection unit 110 randomly selects pixels, but the pixels may not be randomly selected.
  • M pixels may be selected from N pixels in a predetermined area according to a predetermined rule.
  • FIG. 3 is an explanatory diagram showing an example of the processing results of the low-dimensional compression unit 120, the decompression unit 130, and the optimization unit 140 in the image analysis device 10.
  • the expansion unit 130 outputs data such that the phase changes smoothly as compared with the data of the input SAR image (see the central unit in FIG. 3).
  • the optimization unit 140 executes the optimization process using the evaluation function, the observed pixel value is also taken into consideration, so that the processing result reflects the influence of noise to some extent (see the right side in FIG. 3).
  • the spatially correlated phase estimated value converges to the optimum phase such that the phase changes smoothly.
  • This embodiment is suitably applicable to, for example, an application in which a spatially low-frequency phase component is robustly extracted into noise.
  • Spatically low-frequency phase components include, for example, phase delay due to moisture in the atmosphere. Since this phase delay is a phase unrelated to the displacement, the performance of the displacement analysis can be improved by performing estimation according to the present embodiment and excluding it from the observed phase, for example, in the case of displacement analysis.
  • FIG. 4 is a block diagram showing a modified example of the image analysis device of the present embodiment.
  • the image analysis device 10A shown in FIG. 4 further includes a storage unit 131 and a storage unit 141 in addition to the components of the image analysis device 10 shown in FIG.
  • the expansion unit 130 outputs the absolute value and / or the phase of the spatial correlation phase estimation value V samp f s to the storage unit 131.
  • the storage unit 131 stores the absolute value and / or phase of the spatial correlation phase estimation value V samp f s .
  • the expansion unit 130 outputs the absolute value or phase of the spatial correlation phase estimation value V samp f s, or both, to the storage unit 131 each time the loop processing (processing of steps S101 to S104) is executed, and stores the space correlation phase estimation value V samp f s .
  • the unit 131 may store, for example, the absolute value and / or phase of the latest spatially correlated phase estimate V samp f s .
  • the expansion unit 130 may output the absolute value and / or phase of the spatial correlation phase estimation value V samp f s when the loop processing is executed a predetermined number of times.
  • the storage unit 131 stores the data when it is output from the expansion unit 130.
  • the optimization unit 140 outputs the absolute value and / or phase of x samp, s to the storage unit 141.
  • the storage unit 141 stores the absolute value and / or phase of x samp, s .
  • the optimization unit 140 outputs the absolute value and / or phase of x samp, s to the storage unit 141 each time the loop processing is executed, and the storage unit 141 outputs, for example, the latest x samp, s.
  • the absolute value and / or phase of may be preserved in phase.
  • the optimization unit 140 may output the absolute value, phase, or both of x samp, s when the loop processing is executed a predetermined number of times. In that case, the storage unit 141 stores the data when it is output from the optimization unit 140.
  • the absolute value or phase stored in the storage unit 131 or the storage unit 141 can be used for various purposes.
  • the phase is utilized as the phase from which noise has been removed.
  • the absolute value can be used as a value (reliability) indicating the reliability of noise reduction.
  • both the storage unit 131 and the storage unit 141 are provided, only one of the storage unit 131 and the storage unit 141 may be provided. ..
  • FIG. 5 is a block diagram showing a configuration example of the image analysis device of the second embodiment.
  • the image analysis device 20 shown in FIG. 5 further includes a smoothing unit 150 in addition to the components of the image analysis device 10 shown in FIG.
  • the smoothing unit 150 smoothes the output of the low-dimensional compression unit 120.
  • smoothing means smoothing the output of the low-dimensional compression unit 120 in time.
  • steps S101 to S104 and S106 is the same as the processing in the first embodiment.
  • step S104 in the processing of step S104, "f new " is described instead of "fs", but the processing content of step S104 is not changed.
  • the value obtained by smoothing f new using the following equation (2) is defined as fs.
  • f new is the output of the low-dimensional compression unit 120. That is, f new corresponds to f s in the first embodiment.
  • f old is the output of the low-dimensional compression unit 120 in the previous loop processing (processing of steps S101 to S104). That is, f old corresponds to f s calculated by the low-dimensional compression unit 120 in the previous process.
  • f old is stored in the low-dimensional compression unit 120.
  • ⁇ (0 ⁇ ⁇ 1) is the reflection rate of f old .
  • ⁇ is a predetermined value.
  • the pixel selection unit 110 randomly selects pixels.
  • the value of fs calculated in the process of step S104 in one loop process is significantly different from the value of fs calculated in the process of step S104 in the previous loop process. ..
  • FIG. 7 is a block diagram showing a configuration example of the image analysis device of the third embodiment.
  • the configuration of the image analysis device 30 shown in FIG. 7 is basically the same as that of the image analysis device 10 of the first embodiment.
  • the optimization unit 140 is embodied. That is, in the third embodiment, a specific example of the optimization unit 140 is shown.
  • the processing of the pixel selection unit 110, the low-dimensional compression unit 120, and the expansion unit 130 is the same as the processing in the first embodiment.
  • the optimization unit 140 includes a phase averaging unit 142 and an activation unit 143.
  • the phase averaging unit 142 adds the value of the observation signal evaluation function and the value of the spatial correlation evaluation function used in the first embodiment as illustrated in the following equation (3).
  • the second term ( ⁇ means , nm, s) on the right side corresponds to the output of the expansion unit 130. That is, ⁇ corresponds to V samp f s as the spatially correlated phase estimate in the first embodiment.
  • the first term on the right side corresponds to a specific example of the observation signal evaluation function in the first embodiment.
  • nm is synonymous with " samp " in the first embodiment.
  • phase averaging process the process of the phase averaging unit 142 using the data such that the phase output by the expansion unit 130 changes smoothly is referred to as a phase averaging process.
  • the optimization unit 140 initializes the variable i to 0 (step S131). Then, the optimization unit 140 increments the variable i by +1 (step S132). The phase averaging unit 142 in the optimization unit 140 executes the phase averaging process using the equation (3) to obtain a temporary value x temp (step S133).
  • the activation unit 143 performs an activation treatment on x temp to obtain xn m , s.
  • xn m , s is synonymous with x samp, s in the first embodiment.
  • "n" in the xn m , s is represented by the subscript character, but in xn m , s, the "n” in the formula (3) is represented by the subscript character. Same as a variable.
  • step S134 the activation unit 143 non-linearly transforms x temp using the non-linear function g (a) as shown in the following equation (4).
  • I 0 and I 1 are first-order modified Bessel functions of order 0 and order, respectively.
  • the value of g (a) approaches 1 as the value of a increases.
  • the value of g (a) approaches 1 as the value of a becomes larger. That is, in the present embodiment, the maximum value (absolute value) of xn m , s is limited to 1.
  • the activation unit 143 is not limited to the functions exemplified in the equations (5) to (7), and uses another function such that the output value gradually approaches a specific positive value from 0 as the a increases. It is possible.
  • the phase averaging unit 142 and the activation unit 143 may sequentially execute the processes of steps S133 and S134, or simultaneously (parallel) execute the processes of steps S133 and S134 for a plurality of i's. You may.
  • the method of executing in parallel is suitable for parallel calculation by GPU (Graphics Processing Unit), CPU (Central Processing Unit), etc., and can speed up the optimization process.
  • the storage unit 131 and / or the storage unit 141 as shown in FIG. 4 may be added to the image analysis device 30 of the present embodiment.
  • the smoothing unit 150 as shown in FIG. 5 may be added to the image analysis device 30 of the present embodiment.
  • the storage unit 131, the storage unit 141, or both of them, and the smoothing unit 150 may be added to the image analysis device 30 of the present embodiment.
  • the absolute value of the output (xn m , s) of the optimization unit 140 in the present embodiment can be utilized as a value (reliability) indicating the reliability of noise reduction.
  • x This is because it is configured to increase the absolute value of temp . That is, in the equation (3), the value of the first term becomes larger as the phase of the observed pixel value y and the phase of x are closer.
  • the value of the second term is a value obtained by smoothing x calculated in the previous loop processing (processing of steps S101 to S104).
  • the value of x temp becomes large when the value based on the second term and the value based on the first term have a close phase. Therefore, the absolute value of x temp increases as the reliability increases. Therefore, the absolute value of xn m , s also increases as the reliability increases.
  • the absolute values of xn m and s stored in the storage unit 142 can be effectively used as the reliability.
  • FIG. 9 is an explanatory diagram for explaining the effect of the present embodiment.
  • the phase estimation value by the image analysis device 30 is illustrated in the upper right of FIG. Comparing the estimated phase value with the correct answer, there are some points that do not match.
  • the reliability of the phase estimate is illustrated in the lower right of FIG. As mentioned above, the reliability is an absolute value of xn m , s. The higher the absolute value, the higher the reliability.
  • FIG. 10 is an explanatory diagram showing an application example of this embodiment. As an example, it is possible to exclude a region with low reliability in the estimation result and use it for some analysis or the like. In FIG. 10, the black area is an excluded area.
  • the reliability according to this embodiment can be applied to applications such as visually displaying an unreliable region after displacement analysis or the like is executed based on the estimated phase value.
  • FIG. 11 is a block diagram showing a configuration example of the image analysis device of the fourth embodiment.
  • the image analysis device 40 shown in FIG. 11 further includes a front filter unit 160 in addition to the components of the image analysis device 30 shown in FIG. 7.
  • the pixel selection unit 110, the low-dimensional compression unit 120, the expansion unit 130, and the optimization unit 140 including the phase averaging unit 142 and the activation unit 143 are the same as those in the third embodiment. It works in the same way.
  • the phase averaging unit 142 uses an evaluation function different from that of the phase averaging unit 142 in the third embodiment.
  • the pre-filter unit 160 calculates the matrix ⁇ by, for example, the following equation (8).
  • W filter is a predetermined weight.
  • n ⁇ ⁇ (nm) indicates that it is one of the selected M pixels.
  • Equation (8) shows that the numerator is a conjugate complex number of y n, s and y n, s' . Since it contains the product of, it means that the phase difference of the neighboring observation pixels is averaged by the W filter .
  • the phase averaging unit 142 executes the same processing as in step S133 in FIG. However, in the present embodiment, the phase averaging unit 142 calculates x temp using the following equation (9). Other processing is the same as the processing of the image analysis apparatus 30 shown in FIG. Equation (9) corresponds to the equation in which the product of y n , s and the conjugate complex number of y n, s'in equation (3) is replaced by the matrix ⁇ .
  • the front filter unit 160 determines the weight matrix that reflects the average of the phase differences of the neighboring observation pixels, so that the optimum weight can be determined from the information in the observation image. Then, since the phase averaging unit 142 calculates x temp using the evaluation function including the matrix ⁇ reflecting the average of the phase differences of the neighboring observation pixels, the degree of phase noise reduction can be further increased. ..
  • FIG. 12 is a block diagram showing a configuration example of the image analysis device according to the fifth embodiment.
  • the image analysis device 50 shown in FIG. 11 further includes a weight calculation unit 170 in addition to the components of the image analysis device 40 shown in FIG.
  • the weight calculation unit 170 calculates the weight Wn m as shown by the following equation (10). That is, the weight calculation unit 170 calculates the weight Wn m based on the matrix ⁇ .
  • FIG. 13 is a block diagram showing an example of a computer having a CPU.
  • the computer is mounted on the image analysis device of each of the above embodiments.
  • the CPU 1000 realizes each function in the above embodiment by executing a process according to a program (software element: code) stored in the storage device 1001. That is, the pixel selection unit 110, the low-dimensional compression unit 120, in the image analysis devices 10, 10A, 20, 30, 40, 50 shown in FIGS. 1, 4, 5, 7, 11, and 12.
  • the functions of the expansion unit 130, the optimization unit 140, the phase averaging unit 142, the activation unit 143, the smoothing unit 150, the pre-filter unit 160, and the weight calculation unit 170 are realized.
  • a GPU or a combination of a CPU and a GPU can also be used.
  • the storage device 1001 is, for example, a non-transitory computer readable medium.
  • Non-temporary computer-readable media include various types of tangible storage media (tangible storage medium). Specific examples of non-temporary computer-readable media include magnetic recording media (eg, hard disk), optomagnetic recording media (eg, optomagnetic disc), CD-ROM (Compact Disc-Read Only Memory), and CD-R (Compact). Disc-Recordable), CD-R / W (CompactDisc-ReWritable), semiconductor memory (for example, mask ROM, PROM (ProgrammableROM), EPROM (ErasablePROM), flash ROM).
  • the storage device 1001 can also be used as the SAR image storage unit 200 and the weight matrix storage unit 300.
  • the storage device 1001 can also be used as the storage units 131 and 141.
  • the program may be stored in various types of temporary computer readable media (transitory computer readable medium).
  • the temporary computer-readable medium is supplied with a program, for example, via a wired or wireless channel, i.e., via an electrical signal, an optical signal, or an electromagnetic wave.
  • the memory 1002 is realized by, for example, a RAM (RandomAccessMemory), and is a storage means for temporarily storing data when the CPU 1000 executes processing. It is also possible to envision a form in which a program held by the storage device 1001 or a temporary computer-readable medium is transferred to the memory 1002, and the CPU 1000 executes processing based on the program in the memory 1002.
  • a RAM RandomAccessMemory
  • FIG. 14 is a block diagram showing a main part of the image analysis device.
  • the image analysis device 100 shown in FIG. 14 has a pixel selection unit (pixel selection means) 11 (pixel selection in the embodiment) that selects pixels at a plurality of positions in one image of a plurality of images recorded in the same area. (Realized by unit 110) and a low-dimensional compression unit (low-dimensional compression means) 12 (low in the embodiment) that low-dimensionally compresses a complex vector as an evaluation value when the evaluation function is optimized into a low-dimensional space.
  • pixel selection unit pixel selection means
  • low-dimensional compression means 12 low in the embodiment
  • the dimensionality reduction unit 120 returns the compression result of low-dimensional compression unit 12 to the original pixel space, expansion unit for calculating a spatial correlation phase estimate (deployment means) 13 (embodiment, expansion unit
  • the optimization unit (optimization means) that optimizes the evaluation function by bringing the evaluation value closer to the spatial correlation phase estimation value and the pixel at the position selected by the pixel selection unit 11. ) 14 (in the embodiment, it is realized by the optimization unit 140).
  • FIG. 15 is a block diagram showing a main part of an image analysis device of another aspect.
  • the evaluation function used by the optimization unit 14 weights the product of the complex number representing the value of the pixel at the position selected by the pixel selection unit 11 and its conjugate complex number to obtain a plurality of images.
  • the optimization unit 14 includes a term to be added over the above, and the optimization unit 14 performs a non-linear conversion of the value obtained by optimizing the evaluation function and limits the value to a predetermined value or less. It is realized by the conversion unit 143.).
  • Pixel selection means for selecting pixels at a plurality of positions in one image among a plurality of images in which the same area is recorded, and A low-dimensional compression means that low-dimensionally compresses a complex vector as an evaluation value when the evaluation function is optimized into a low-dimensional space, A decompression means for returning the compression result by the low-dimensional compression means to the original pixel space and calculating a spatial correlation phase estimate, and a decompression means.
  • An image analysis apparatus including an optimization means for optimizing the evaluation function by bringing the evaluation value closer to the spatial correlation phase estimation value and the pixel at a position selected by the pixel selection means.
  • the evaluation function includes a term that weights the product of a complex number representing the value of a pixel at a position selected by the pixel selection means and its conjugate complex number and adds them over the plurality of images.
  • the image analysis apparatus according to Appendix 1, wherein the optimization means includes a limiting means that non-linearly transforms a value obtained by optimizing the evaluation function and limits it within a predetermined value.
  • the low-dimensional compression means multiplies the complex vector by a matrix for low-dimensional compression.
  • the expansion means is an image analysis device according to Appendix 1 or Appendix 2 for deriving the spatial correlation phase estimate by multiplying the compression result by the low-dimensional compression means by a matrix for expansion.
  • the supplementary note further includes a smoothing means (in the embodiment, realized by the smoothing unit 150) that smoothes the compression result by the low-dimensional compression means with the past compression result and then supplies the compression result to the decompression means.
  • a smoothing means in the embodiment, realized by the smoothing unit 150
  • An image analysis device according to any one of 1 to 3.
  • Supplementary Note 5 An image analysis apparatus according to any one of Supplementary note 1 to Supplementary note 4, further comprising an absolute value storage means (in the embodiment, realized by the storage unit 131) for storing the absolute value of the output of the expansion means. ..
  • Supplementary note 6 An image analysis apparatus according to any one of Supplementary note 1 to Supplementary note 5, further comprising a phase storage means (in the embodiment, realized by the storage unit 131) for storing the phase of the output of the expansion means.
  • a phase storage means in the embodiment, realized by the storage unit 131) for storing the phase of the output of the expansion means.
  • Supplementary note 7 Image analysis of any one of Supplementary note 1 to Supplementary note 4 further comprising an absolute value storage means (in the embodiment, realized by the storage unit 141) for storing the absolute value of the output of the optimization means. Device.
  • Supplementary note 8 An image of any one of Supplementary note 1 to Supplementary note 4 or Supplementary note 7 further comprising a phase saving means (in the embodiment, realized by the storage unit 141) for storing the phase of the output of the optimization means. Analytical device.
  • the evaluation function includes a term added over the plurality of images by weighting the product of a complex number representing the value of the pixel at a position selected by the pixel selection means and the conjugate complex number thereof.
  • the image analysis method of Appendix 9 in which the value obtained by optimizing the evaluation function is non-linearly converted and limited to a predetermined value.
  • the evaluation function includes a term added over the plurality of images by weighting the product of the complex number representing the value of the pixel at the position selected by the pixel selection means and the conjugate complex number thereof.
  • the computer-readable recording medium of Appendix 11 containing an image analysis program that causes a computer to perform a process of non-linearly converting the value obtained by optimizing the evaluation function and limiting the value within a predetermined value.
  • the evaluation function includes a term added over the plurality of images by weighting the product of a complex number representing the value of the pixel at a position selected by the pixel selection means and the conjugate complex number thereof.
  • the image analysis program of Appendix 13 for causing a computer to perform a process of non-linearly converting the value obtained by optimizing the evaluation function and limiting it within a predetermined value.
  • Pixel selection unit 12 Low-dimensional compression unit 13
  • Expansion unit 14 Optimization unit 15 Restriction unit 10, 10A, 20, 30, 40, 50 Image analysis device 100, 101 Image analysis device 110 Pixel selection unit 120 Low-dimensional compression unit 131, 141 Storage unit 130 Expansion unit 140 Optimization unit 142 Phase averaging unit 143 Activation unit 150 Smoothing unit 160
  • Front filter unit 170 Weight calculation unit 200
  • SAR image storage unit 300 Weight matrix storage unit 1000

Abstract

画像解析装置１００は、同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する画素選択部１１と、評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する低次元圧縮部１２と、低次元圧縮部１２による圧縮結果を元の画素空間に戻し、空間相関位相推定値を算出する展開部１３と、評価値を、空間相関位相推定値と、画素選択部１１によって選択された位置の画素とに近づけることによって、評価関数を最適化する最適化部１４とを含む。

Description

画像解析装置および画像解析方法
　本発明は、合成開口レーダの受信電磁波の位相に基づいて画像解析を行う画像解析システムに適用可能な画像解析装置および画像解析方法に関する。
　合成開口レーダ（ＳＡＲ：Synthetic Aperture Radar）技術は、人工衛星や航空機などの飛翔体が移動しながら電磁波を送受信し、大きな開口を持つアンテナと等価な観測画像を得るための技術である。合成開口レーダは、例えば、地表からの反射波を信号処理して、標高や地表変位等を解析するために利用される。特に、精度が必要である場合等に、解析装置は、合成開口レーダによって得られる時系列のＳＡＲ画像（ＳＡＲデータ）を入力とし、入力されたＳＡＲ画像を時系列解析する。
　標高や地表変位を解析するための有効な手法として、干渉ＳＡＲ解析がある。干渉ＳＡＲ解析では、違う時期に撮影された複数（例えば、２枚）のＳＡＲ画像を構成する電波信号の位相差が計算される。そして、位相差に基づいて撮影時期間で生じた飛翔体と地面間の距離の変化が検出される。
　特許文献１には、コヒーレンス行列を用いる解析手法が記載されている。コヒーレンス行列は、複数の複素画像の同じ位置にあたる画素の相関を表す。
　コヒーレンスは、Ｓ（Ｓ≧２）枚のＳＡＲ画像のうちの複数のＳＡＲ画像における同じ位置にあたる画素の複素相関で計算される。ＳＡＲ画像のペアを（ｐ、ｑ）とし、コヒーレンス行列の成分をｃｐ，ｑとする。ｐ，ｑは、それぞれ、Ｓ以下の値であり、Ｓ枚のＳＡＲ画像のいずれかを示す。ＳＡＲ画像のペアについて、位相θｐ，ｑ（具体的には、位相差）が算出される。そして、コヒーレンス算出対象の画素を含む所定領域内の複数の画素についてｅｘｐ（－ｊθｐ，ｑ）が平均化された値の絶対値が、コヒーレンス行列の成分ｃｐ，ｑとなる。
　ｃｐ，ｑの絶対値｜ｃｐ，ｑ｜から、位相θｐ，ｑの分散の大きさを把握可能である。
　コヒーレンス行列は、分散のように位相ノイズ量の程度を推測可能な情報を含む。
　地表等の変位解析のために、位相θｐ，ｑが変位速度および撮影時刻差に相関することが利用される。例えば、位相差の平均値に基づいて変位が推定される。なお、位相ノイズ量を利用して変位解析の精度を検証することが可能である。したがって、コヒーレンス行列は、変位解析のために使用可能である。
　標高解析のために、位相θｐ，ｑが解析対象物の標高および飛翔体間距離（例えば、飛翔体の２つの撮影位置の間の距離）に相関することが利用される。例えば、位相差の平均値に基づいて標高が推定される。なお、位相ノイズ量を利用して標高解析の精度を検証することが可能である。したがって、コヒーレンス行列は、標高解析のために使用可能である。
　特許文献１には、コヒーレンス行列に対して変位等のモデルをフィッティングし、ノイズの影響を除外した位相を復元する手法が記載されている。また、非特許文献１にも、類似の手法が開示されている。
国際公開第２０１１／００３８３６号
A. M. Guarnieri et.al., "On the Exploitation of Target Statistics for SAR Interferometry Applications", IEEE Transactions on Geoscience and Remote Sensing, Vol. 46, No. 11, pp.3436-3443, November 2008
　特許文献１に記載された手法および非特許文献１に記載された手法によれば、位相差に含まれるノイズを減らすことができる。しかし、変位解析および標高解析が実行されるときに、位相差に含まれるノイズはできるだけ少ないことが好ましいので、より多くのノイズが除去されることが要望される。
　本発明は、位相ノイズ低減の程度をより大きくすることができる画像解析装置および画像解析方法を提供することを目的とする。
　本発明による画像解析装置は、同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する画素選択手段と、評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する低次元圧縮手段と、低次元圧縮手段による圧縮結果を元の画素空間に戻し、空間相関位相推定値を算出する展開手段と、評価値を、空間相関位相推定値と、画素選択手段によって選択された位置の画素とに近づけることによって、評価関数を最適化する最適化手段とを含む。
　本発明による画像解析方法は、同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択し、評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮し、低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出し、評価値を、空間相関位相推定値と、選択された位置の画素とに近づけることによって、評価関数を最適化する。
　本発明による画像解析プログラムは、コンピュータに、同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する処理と、評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する処理と、低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出する処理と、評価値を、空間相関位相推定値と、選択された位置の画素とに近づけることによって、評価関数を最適化する処理とを実行させる。
　本発明によれば、位相ノイズ低減の程度をより大きくすることができる。
第１の実施形態の画像解析装置を示すブロック図である。 第１の実施形態の画像解析装置の動作を示すフローチャートである。 低次元圧縮部、展開部、および最適化部の処理結果の一例を示す説明図である。 第１の実施形態の画像解析装置の変形例を示すブロック図である。 第２の実施形態の画像解析装置を示すブロック図である。 第２の実施形態の画像解析装置の動作を示すフローチャートである。 第３の実施形態の画像解析装置の構成例を示すブロック図である。 第３の実施形態の画像解析装置における最適化部の動作を示すフローチャートである。 第３の実施形態の効果を説明するための説明図である。 第３の実施形態の応用例を示す説明図である。 第４の実施形態の画像解析装置の構成例を示すブロック図である。 第５の実施形態の画像解析装置の構成例を示すブロック図である。 ＣＰＵを有するコンピュータの一例を示すブロック図である。 画像解析装置の主要部を示すブロック図である。 他の態様の画像解析装置の主要部を示すブロック図である。
　以下、本発明の実施形態を図面を参照して説明する。
　下記の実施形態では、画像解析装置は、重みづけられた観測信号（観測画素値）とノイズがない画素値とに基づく第１の評価関数（観測信号評価関数）と、画像内で位相差が滑らかに分布することを表す第２の評価関数（空間相関評価関数）とを使用する。そして、画像解析装置は、観測信号評価関数と空間相関評価関数とを統合した評価関数を使用して、ノイズがないと推定される位相を算出する。「統合」は、例えば「和」である。以下、単に「評価関数」と表現された場合に、それは、観測信号評価関数と空間相関評価関数とが統合された評価関数を意味する。
　本実施形態の画像解析装置は、重みとしてコヒーレンス行列の絶対値の逆行列に負号を付けたものを用いる。ただし、コヒーレンス行列の絶対値の逆行列に負号を付けたものは一例であり、画像解析装置は、重みとして、他のパラメータを用いてもよい。重みとしてコヒーレンス行列の絶対値の逆行列に負号を付けたものが用いられる場合に、画像解析装置は、コヒーレンス行列算出部を含んでいてもよい。コヒーレンス行列算出部は、例えば、以下のようにコヒーレンス行列を算出する。
　すなわち、コヒーレンス行列算出部は、例えばＳＡＲ画像格納部に格納されているＳ枚のＳＡＲ画像（複素画像：振幅および位相情報を含む。）を対象としてコヒーレンス行列Ｃを算出する。例えば、ＳＡＲ画像のペアを（ｐ，ｑ）とし、コヒーレンス行列Ｃの成分をｃｐ，ｑとする。ｐ，ｑは、それぞれ、Ｓ以下の値であり、Ｓ枚のＳＡＲ画像のいずれかを示す。コヒーレンス行列算出部は、ＳＡＲ画像のペアについて、位相θｐ，ｑ（具体的には、位相差）を算出する。そして、コヒーレンス行列算出部は、コヒーレンス算出対象の画素を含む所定領域内の複数の画素についてｅｘｐ（－ｊθｐ，ｑ）を平均化した値を、コヒーレンス行列Ｃの成分ｃｐ，ｑとする。
　また、ＳＡＲ画像ｐにおける強度をＡp、ＳＡＲ画像ｑにおける強度をＡqとして、Ａp・Ａq・ｅｘｐ（－ｊθp, q）を平均化してもよい。また、Ａp・Ａq・ｅｘｐ（－ｊθp, q）の平均として得られた行列に対して、その対角成分の総和をＳで割った値によって、行列の要素それぞれを割ってもよい。また、Ａp・Ａq・ｅｘｐ（－ｊθp, q）の平均として得られた行列に対して、その対角成分の-1/2乗を対角成分として持つ対角行列を左右から乗算してもよい。
　なお、このようなコヒーレンス行列の算出方法は一例であり、コヒーレンス行列算出部は、種々の手法でコヒーレンス行列を算出できる。
　観測信号評価関数として、様々に表現可能な関数のいずれも採用しうる。以下の実施形態では、観測信号評価関数は、一例として、少なくとも、評価値としての、ノイズがない画素値ｘすなわち推定される画素値ｘに加えて、重みＷ（以下の実施形態では、コヒーレンス行列Ｃの絶対値|Ｃ|の逆行列|Ｃ|-1に負号を付けたものである－|Ｃ|-1）、および観測画素値ｙをパラメータとする評価関数（例えば、積和演算式）である。
　画像解析装置は、空間相関予測部を含んでいてもよい。空間相関予測部は、画像における画素の相関係数を予測し、相関係数を要素とする行列である空間相関行列Ｋを生成する。空間相関行列Ｋは、画素の相関係数を要素とする行列である。ただし、空間相関行列Ｋは空間相関の表現の一例であって、他の表現が用いられてもよい。
　例えば、空間相関予測部は、ＳＡＲ画像における解析対象領域（全領域でもよい。）に関する事前情報（既知のデータ）に基づいて、空間相関（解析対象領域における画素の相関）を作成（算出）する。事前情報は、事前に空間相関予測部に与えられる。事前情報として、例えば、以下のような情報が例示される。
・画像解析装置１０が適用された画像解析システムの解析対象物（例えば、変位解析の対象である構造物）の概略形状。一例として、鉄塔などの細長い形状の構造物が対象物である場合には、構造物の延伸方向に画素が滑らかな相関を持つように空間相関を作成する。
・相関が既知であるグラフ（ノード：画素、エッジ（重み）：相関）における重み。この事前情報を使用する場合、重みに基づく値を要素とする空間相関行列Ｋを作成する。
・過去に得られたＳＡＲ画像から得られる情報。すなわち、解析対象領域の既知画像における画素の相関。
　事前情報に基づいて生成される空間相関行列Ｋを含む空間相関評価関数は、画像内で位相差が滑らかに分布することを表す関数に相当する。
　なお、事前情報は、上記の例示に限られず、解析対象物の空間相関を予測可能であれば、他の類いの情報でもよい。
　空間相関評価関数として、様々に表現可能な関数のいずれも採用しうる。以下の実施形態では、空間相関評価関数は、一例として、少なくとも、評価値としての画素値ｘに加えて、一種の重みに相当する空間相関行列Ｋ（実際には、空間相関行列Ｋの逆行列Ｋ-1）をパラメータとする評価関数（例えば、積和演算式）である。
　そして、以下の実施形態では、画像解析装置は、所定領域内の全画素（画素１～画素Ｎ）について観測信号評価価数を用いた最適化処理を実行するのではなく、全画素から選択されたＭ個（Ｍ＜Ｎ）画素について観測信号評価価数を用いた最適化処理を実行する。
　また、上述したように評価関数には観測信号評価関数と空間相関評価関数とが含まれるが、画像解析装置は、空間相関評価関数における空間相関行列Ｋを低次元分解（低ランク近似）して最適化処理を実行する。なお、詳しくは、画像解析装置は、空間相関行列Ｋの逆行列Ｋ-1を実質的に低次元分解する。画像解析装置は、例えば、Nystrom 近似と呼ばれる手法を用いて低次元分解を行うことができる。その場合、画像解析装置は、（１）式を用いて空間相関行列Ｋの逆行列Ｋ-1を低次元分解する。
Figure JPOXMLDOC01-appb-M000001
　（１）式において、Ｇは、ｄ×ｄの行列である（ｄ＜Ｎ）。Ｖは、Ｓ×ｄの行列（Ｖ＝ｖ１，ｖ２，・・・，ｖｄ）である。Λは、Ｎ×Ｎの対角行列である。
実施形態１．
　図１は、第１の実施形態の画像解析装置の構成例を示すブロック図である。図１に示す画像解析装置１０は、画素選択部１１０、低次元圧縮部１２０、展開部１３０、および最適化部１４０を含む。なお、上述したように、画像解析装置１０が、コヒーレンス行列算出部および空間相関予測部を含むことがある（図１において図示せず）。
　ＳＡＲ画像格納部２００には、あらかじめＳＡＲ画像が格納される。重み行列格納部３００には、最適化処理で用いられる重みがあらかじめ格納される。本実施形態では、重み行列格納部３００に、重みとしての重み行列が格納される。なお、ＳＡＲ画像格納部２００および重み行列格納部３００は、画像解析装置１０の外部に設置されていてもよいし、画像解析装置１０に含まれていてもよい。
　画素選択部１１０は、同一地域が記録された複数のＳＡＲ画像から、ある１枚のＳＡＲ画像中の画素を複数選択する。そして、画素選択部１１０は、他のＳＡＲ画像における、上記の１枚のＳＡＲ画像中の画素と同一の位置の画素を特定する。画素選択部１１０は、特定した画素の位置を指定する。低次元圧縮部１２０は、画素選択部１１０が指定した位置の画素を低次元空間に写像する。展開部１３０は、画素選択部１１０が指定した位置の画素を低次元空間から元の画素空間に写像する。
　最適化部１４０は、上記の１枚のＳＡＲ画像と他のＳＡＲ画像とのペアについての重み行列を用いる。なお、画素選択部１１０は、１枚のＳＡＲ画像と他のＳＡＲ画像との中での共通した位置を指定する。そして、最適化部１４０は、展開部１３０から出力された演算結果と画素選択部１１０が指定した位置の観測画素の画素値（観測画素値）との両方に対して近くなるような位相を持つ画素（複素ベクトルで表される。後述するｘsamp,sに相当する。）を導出する。
　なお、最適化部１４０は、空間相関評価関数を最適化（例えば、最大化）することによって、位相を持つ複素ベクトルを展開部１３０から出力された演算結果に近づける。また、最適化部１４０は、観測信号評価関数を最適化（例えば、最大化）することによって、複素ベクトルを画素選択部１１０が指定した位置の観測画素の画素値に近づける。また、観測信号評価関数が最適化されるときに、複素ベクトルは、ノイズが比較的小さい領域における観測画素の画素値に近づく。
　また、本実施形態では、画素選択部１１０は、画素の位置を指定したが、ＳＡＲ画像格納部２００から該当する画素（具体的には、画素値）を読み出して、読み出した画素を、低次元圧縮部１２０、展開部１３０、および最適化部１４０に供給してもよい。
　次に、図２のフローチャートを参照して本実施形態の画像解析装置１０の動作を説明する。なお、重み行列格納部３００には、重みＷがあらかじめ格納される。また、（１）式で表される空間相関行列Ｋの逆行列Ｋ-1が画像解析装置１０にあらかじめ設定される。なお、低次元圧縮部１２０は、空間相関行列Ｋの逆行列Ｋ-1の一部である行列ＧＶTを使用する。また、展開部１３０は、空間相関行列Ｋの逆行列Ｋ-1の一部である行列Ｖを使用する。よって、行列ＧＶTと行列Ｖとが、画像解析装置１０に設定されてもよい。
　画素選択部１１０は、Ｍ個（Ｍ＜Ｎ）の乱数を生成する（ステップＳ１０１）。Ｎは、ＳＡＲ画像の総画素数を示す。なお、Ｍは２以上の値である。
　以下、選択される画素の画素番号をｎｍとする（ｍ：１～Ｍ）。なお、ｎｍは、１～Ｍのうちのいずれかになる。Ｓは、総ＳＡＲ画像数を示す。画像番号をｓとする（ｓ：１～Ｓ）。また、ｎｍをｓａｍｐと表現することがある。画素選択部１１０は、生成した複数の乱数のそれぞれを、選択される画素の画素番号とする。すなわち、画素選択部１１０は、ランダムに画素を選択する。
　展開部１３０は、前回のｆs にＶsampを乗算する（ステップＳ１０２）。前回のｆs は、直前に実行されたステップＳ１０４の処理の演算結果を意味する。Ｖsampのそれぞれ（samp：１～Ｍ）は、行列Ｖから、ｎ１行、ｎ２行、・・・、ｎＭ行を取り出して生成される行ベクトルである。すなわち、展開部１３０は、画素番号ｎ１，ｎ２，・・・，ｎＭに対応する行のそれぞれの画素について、行列Ｖにおける対応する行ベクトルを用いてＶsampｆs を算出する。以下、Ｖsampｆs を、空間相関位相推定値ということがある。
　なお、後述するように、ｆs は、上記の（１）式におけるＧＶTにｘsamp,sが乗算されて得られた値（複素数のデータ）に相当する。よって、Ｖsampｆs は、（１）式におけるＶＧＶTにｘsamp,sが乗算されて得られた値（複素数のデータ）に相当する。"samp,s"における"s"は、画像番号であるから、例えば、"samp,1"は、画像番号１のＳＡＲ画像における選択された画素を指定するインデックスである。
　最適化部１４０は、評価関数を最適化（例えば、最大化）する（ステップＳ１０３）。評価関数が最適化されることによって、ｘsamp,sは最適値に近づく。上述したように、評価関数は、少なくとも評価値としてのノイズがないと想定された画素値ｘに加えて、重みＷ、および観測画素値ｙをパラメータとして含む観測信号評価関数と、空間相関行列Ｋ（実際には、空間相関行列Ｋの逆行列Ｋ-1）をパラメータとする空間相関評価関数とが統合された評価関数である。
　最適化部１４０は、ステップＳ１０３の処理で、ｘsamp,sを、観測信号評価関数に関して、Ｗsamp（本実施形態では、選択された画素に対応するＳ×Ｓ個のコヒーレンス行列の逆行列）で重み付けられたｙsamp,sに近づけ、かつ、空間相関評価関数に関して、空間相関位相推定値であるＶsampｆs に近づける処理を実行することになる。
　低次元圧縮部１２０は、最適化部１４０によって最適化されたｘsamp,sを、ＧＶT sampを用いて圧縮する（ステップＳ１０４）。すなわち、低次元圧縮部１２０は、ｘsamp,sとＧＶT sampとを乗算する。低次元圧縮部１２０は、さらに、乗算結果を（Ｎ／Ｍ）倍する。乗算結果が（Ｎ／Ｍ）倍された値を、ｆs と表現する。なお、Ｎ個の全画素から選択されたＭ個の画素を対象としてステップＳ１０２～Ｓ１０４の処理が実行されるので、ｘsamp,sとＧＶT sampとの乗算結果は、本来の値の（Ｍ／Ｎ）倍になっている。低次元圧縮部１２０は、乗算結果を（Ｎ／Ｍ）倍して本来の値に戻す。
　例えば、低次元圧縮部１２０は、メモリを内蔵する。低次元圧縮部１２０は、ｆsをメモリ（図１において図示せず）に一時保存して展開部１３０に渡すように構成されている。
　終了条件が成立した場合には、処理を終了する（ステップＳ１０６）。終了条件が成立していない場合には、ステップＳ１０１に戻る。なお、画素選択部１１０は、次に実行されるステップＳ１０１の処理において、前回に選択した画素群とは別の画素群（一部重複していてもよい。）を選択する。
　終了条件は、例えば、展開部１３０から出力される空間相関位相推定値としてのＶsampｆs が、前回のステップＳ１０２の処理で得られた値に対して、あらかじめ定められた所定値未満であることである。すなわち、終了条件は、Ｖsampｆs の変動量が小さくなっていることである。換言すれば、終了条件は、Ｖsampｆs が最適値に収束したと判定されたことである。また、終了条件は、最適化部１４０から出力されるｘsamp,sが、前回のステップＳ１０３の処理で得られた値に対して、あらかじめ定められた所定値未満であることでもよい。すなわち、終了条件は、ｘsamp,sの変動量が小さくなっていることでもよい。換言すれば、終了条件は、ｘsamp,sが最適値に収束したと判定されたことでもよい。なお、そのような終了条件は一例であって、他の条件、例えば、Ｓ１０１～Ｓ１０４の処理が所定回実行されるという条件が用いられてもよい。
　以上に説明したように、本実施形態では、低次元圧縮部１２０は、ＧＶT sampｘsamp,sを算出する。展開部１３０は、ＧＶT sampｘsamp,sが（Ｎ／Ｍ）倍された値（ｆs ）をＶによって展開する。最適化部１４０は、観測信号評価関数に関して、前回選択された画素群とは異なる画素群を用いて空間相関評価関数を最適化する。低次元圧縮部１２０、展開部１３０および最適化部１４０が処理を繰り返すことによって、ｘsamp,sは最適値に収束する。
　すなわち、画像解析装置１０は、空間相関評価関数を最適化することによって位相差が滑らかに分布し、かつ、観測信号評価関数を最適化することによって観測信号の位相に近い位相を得る。よって、位相ノイズ低減の程度をより大きくすることができる。
　また、空間相関行列Ｋの逆行列Ｋ-1が低次元分解されているので、低次元分解されていない場合に比べてメモリの使用量を削減できる。具体的には、メモリの使用量は、（ｆs の次元数×画像数）＋選択された画素の数×画像数２）である。
　また、最適化部１４０は、観測信号評価関数に関して、全画素からランダムに選択された画素を対象として最適化処理を実行する。よって、最適化部１４０の演算量が削減される。
　なお、本実施形態では、画素選択部１１０は、ランダムに画素を選択するが、ランダムに選択しなくてもよい。例えば、あらかじめ定められている所定のルールに従って、所定の領域内のＮ個の画素からＭ個の画素を選択してもよい。
　図３は、画像解析装置１０における低次元圧縮部１２０、展開部１３０、および最適化部１４０の処理結果の一例を示す説明図である。図３の左側には、ノイズが含まれる３枚の入力ＳＡＲ画像が例示されている。例えば、縦方向に帯状に伸びる領域が、比較的高いノイズが含まれる領域である。展開部１３０は、入力ＳＡＲ画像のデータと比較して、位相が滑らかに変化するようなデータを出力する（図３における中央部参照）。最適化部１４０が評価関数を用いて最適化処理を実行すると、観測画素値も考慮されるので、処理結果には、ある程度のノイズの影響が反映される（図３における右側参照）。
　しかし、低次元圧縮部１２０および展開部１３０による処理が繰り返されると、空間相関位相推定値は、位相が滑らかに変化するような最適な位相に収束する。
　本実施形態は、例えば、空間的に低周波の位相成分をノイズに頑健に抽出するような用途に好適に適用可能である。空間的に低周波の位相成分として、例えば、大気中の水分による位相遅延がある。この位相遅延は、変位と無関係な位相であるので、例えば変位解析の際には本実施形態による推定を行い、観測された位相から除外することによって変位解析の性能を向上させることができる。
［変形例］
　図４は、本実施形態の画像解析装置の変形例を示すブロック図である。図４に示す画像解析装置１０Ａは、図１に示された画像解析装置１０の構成要素に加えて、保存部１３１と保存部１４１とをさらに含む。
　画像解析装置１０Ａにおいて、展開部１３０は、空間相関位相推定値Ｖsampｆs の絶対値もしくは位相またはそれら双方を保存部１３１に出力する。保存部１３１は、空間相関位相推定値Ｖsampｆs の絶対値もしくは位相またはそれら双方を保存する。展開部１３０が、１回のループ処理（ステップＳ１０１～Ｓ１０４の処理）が実行される度に空間相関位相推定値Ｖsampｆs の絶対値もしくは位相またはそれら双方を保存部１３１に出力し、保存部１３１が、例えば、最新の空間相関位相推定値Ｖsampｆs の絶対値もしくは位相またはそれら双方を保存するようにしてもよい。また、展開部１３０が、所定回のループ処理が実行された場合に、空間相関位相推定値Ｖsampｆs の絶対値もしくは位相またはそれら双方を出力するようにしてもよい。その場合、保存部１３１は、展開部１３０からデータが出力されたときにそれを保存する。
　また、最適化部１４０は、ｘsamp,sの絶対値もしくは位相またはそれら双方を保存部１４１に出力する。保存部１４１は、ｘsamp,sの絶対値もしくは位相またはそれら双方を保存する。最適化部１４０が、１回のループ処理が実行される度にｘsamp,sの絶対値もしくは位相またはそれら双方を保存部１４１に出力し、保存部１４１が、例えば、最新のｘsamp,sの絶対値もしくは位相または双方を位相を保存するようにしてもよい。また、最適化部１４０が、所定回のループ処理が実行された場合に、ｘsamp,sの絶対値もしくは位相または双方を出力するようにしてもよい。その場合、保存部１４１は、最適化部１４０からデータが出力されたときにそれを保存する。
　保存部１３１または保存部１４１に保存された絶対値または位相は、種々の用途に使用可能である。例えば、位相は、ノイズが除去された位相として活用される。また、絶対値を、ノイズ除去の信頼性を表す値（信頼度）として使用可能である。
　なお、図４には、保存部１３１と保存部１４１との双方が設けられていることが例示されているが、保存部１３１と保存部１４１とのうちの一方のみが設けられていてもよい。
実施形態２．
　図５は、第２の実施形態の画像解析装置の構成例を示すブロック図である。図５に示す画像解析装置２０は、図１に示された画像解析装置１０の構成要素に加えて、平滑化部１５０をさらに含む。平滑化部１５０は、低次元圧縮部１２０の出力を平滑化する。なお、本実施形態において、平滑化は、低次元圧縮部１２０の出力を時間的に滑らかにすることを意味する。
　次に、図６のフローチャートを参照して本実施形態の画像解析装置２０の動作を説明する。ステップＳ１０１～Ｓ１０４，Ｓ１０６の処理は、第１の実施形態における処理と同じである。図６では、ステップＳ１０４の処理において、"ｆs "ではなく"ｆnew"と表記されているが、ステップＳ１０４の処理内容が変更されているわけではない。
　平滑化部１５０は、例えば下記の（２）式を用いて、ｆnewを平滑化した値をｆs とする。（２）式において、ｆnewは、低次元圧縮部１２０の出力である。すなわち、ｆnewは、第１の実施形態におけるｆs に相当する。ｆoldは、１回前のループ処理（ステップＳ１０１～Ｓ１０４の処理）での低次元圧縮部１２０の出力である。すなわち、ｆoldは、前回の処理で低次元圧縮部１２０が算出したｆs に相当する。ｆoldは、低次元圧縮部１２０において保存されている。ρ（０＜ρ＜１）は、ｆoldの反映率である。ρは、あらかじめ定められた値である。
　ｆs ＝（１－ρ）ｆold＋ρｆnew　　　　　　　　（２）
　第１の実施形態では、各々のループ処理において、画素選択部１１０がランダムに画素を選択する。その結果、１回のループ処理におけるステップＳ１０４の処理で算出されたｆs の値が、前回のループ処理におけるステップＳ１０４の処理で算出されたｆs の値から大幅に変わっている可能性がある。
　本実施形態では、展開部１３０に渡されるｆs の値が前回渡されたｆs の値に近くなることが期待される。よって、早めに適切な空間相関位相推定値に収束する可能性が高まる。
実施形態３．
　図７は、第３の実施形態の画像解析装置の構成例を示すブロック図である。図７に示す画像解析装置３０の構成は、基本的に、第１の実施形態の画像解析装置１０と同じである。ただし、第３の実施形態では、最適化部１４０が具体化されている。すなわち、第３の実施形態では、最適化部１４０の一具体例が示される。画素選択部１１０、低次元圧縮部１２０および展開部１３０の処理は、第１の実施形態における処理と同じである。
　図７に示すように、本実施形態では、最適化部１４０は、位相平均化部１４２と活性化部１４３とを含む。位相平均化部１４２は、第１の実施形態で使用された観測信号評価関数の値と空間相関評価関数の値とを、下記の（３）式に例示するように加算する。
Figure JPOXMLDOC01-appb-M000002
　（３）式において、右辺の第２項（χmeans,nm,s）は、展開部１３０の出力に相当する。すなわち、χは、第１の実施形態における空間相関位相推定値としてのＶsampｆs に相当する。
　（３）式において、右辺の第１項は、第１の実施形態における観測信号評価関数の一具体例に相当する。第１項において、ｓおよびｓ’は、Ｓ枚のＳＡＲ画像から選択される画像の画像番号に相当する。なお、上述したように、ｎｍは、第１の実施形態における"ｓａｍｐ"と同義である。
　次に、図８のフローチャートを参照して本実施形態における最適化部１４０の動作を説明する。以下、展開部１３０が出力する位相が滑らかに変化するようなデータを用いる位相平均化部１４２の処理を、位相平均化処理と呼ぶ。
　最適化部１４０は、変数ｉを０に初期化する（ステップＳ１３１）。そして、最適化部１４０は、変数ｉを＋１する（ステップＳ１３２）。最適化部１４０における位相平均化部１４２は、（３）式を用いた位相平均化処理を実行し、一時的な値であるｘtempを得る（ステップＳ１３３）。
　活性化部１４３は、ｘtempに活性化処理を施して、ｘnm,sを得る。（ステップＳ１３４）。なお、ｘnm,sは、第１の実施形態におけるｘsamp,sと同義である。また、（３）式では、ｘnm,sにおける"n"が下付文字で表されているが、ｘnm,sは、（３）式における"n"が下付文字で表されている変数と同じである。
　活性化部１４３は、ステップＳ１３４の処理で、下記の（４）式のように、非線形関数ｇ（ａ）を用いて、ｘtempを非線形変換する。
Figure JPOXMLDOC01-appb-M000003
　非線形関数ｇ（ａ）として、以下のような関数を使用可能である。
　ｇ（ａ）＝Ｉ1（２ａ）／Ｉ0（２ａ）　　　　　　　（５）
Figure JPOXMLDOC01-appb-M000004
　ｇ（ａ）＝ｔａｎｈ（ａ）　　　　　　　　　　　　（７）
　（５）式におけるＩ０およびＩ１は、それぞれ、０次と１次の第一種変形ベッセル関数である。（５）式を用いる場合、ｇ（ａ）の値は、ａの値が大きくなるほど１に近づく。（６）式を用いる場合にも（７）式を用いる場合にも、ｇ（ａ）の値は、ａの値が大きくなるほど１に近づく。つまり、本実施形態では、ｘnm,sの最大値（絶対値で）は、１に制限される。
　なお、（５）～（７）式に例示した関数に限られず、活性化部１４３は、出力値が、ａの増加に伴って０から特定の正値に漸近するような他の関数を使用可能である。
　ステップＳ１３４の処理が実行された後、最適化部１４０は、ｉ＝Ｓになっているか否か確認する（ステップＳ１３５）。すなわち、全ＳＡＲ画像について、ステップＳ１３３，Ｓ１３４の処理が実行されたか否か確認する。ｉ＝Ｓになっている場合には、処理を終了する。そうでない場合には、ステップＳ１３２に戻る。
　なお、位相平均化部１４２および活性化部１４３が、逐次的にステップＳ１３３，Ｓ１３４の処理を実行してもよいし、複数のｉについて同時に（並列的に）ステップＳ１３３，Ｓ１３４の処理を実行してもよい。並列的に実行する方式は、ＧＰＵ（Graphics Processing Unit）やＣＰＵ（Central Processing Unit）等による並列計算に適し、最適化処理を高速化することができる。
　また、本実施形態の画像解析装置３０に、図４に示されたような保存部１３１もしくは保存部１４１またはそれら双方が追加されてもよい。また、本実施形態の画像解析装置３０に、図５に示されたような平滑化部１５０が追加されてもよい。また、本実施形態の画像解析装置３０に、保存部１３１もしくは保存部１４１またはそれら双方と平滑化部１５０とが追加されてもよい。
　本実施形態における最適化部１４０の出力（ｘnm,s）の絶対値を、ノイズ除去の信頼性を表す値（信頼度）として活用可能である。画像解析装置３０が、空間的滑らかさの程度と、評価値としてのｘ（ｘnm,s）と観測画素値ｙ（ｙnm,s）との合致度との双方が高くなる場合に、ｘtempの絶対値を大きくするように構成されているからである。すなわち、（３）式において、第１項の値は、観測画素値ｙの位相とｘの位相とが近いほど大きくなる。また、第２項の値は、前回のループ処理（ステップＳ１０１～Ｓ１０４の処理）で計算されたｘが滑らかにされた値である。第２項の値と第１項の値とが加算されると、第２項に基づく値と第１項に基づく値とが近い位相を持つ場合に、ｘtempの値は大きくなる。よって、ｘtempの絶対値は、信頼度が高いほど大きくなる。よって、ｘnm,sの絶対値も、信頼度が高いほど大きくなる。
　なお、特に、保存部１４２が設けられている場合には、保存部１４２に保存されたｘnm,sの絶対値を信頼度として有効活用できる。
　図９は、本実施形態の効果を説明するための説明図である。画像解析装置３０による位相の推定値が、図９における右上に例示されている。位相の推定値と正解とを比較すると、一致していない箇所が存在する。位相の推定値の信頼度が、図９における右下に例示されている。上述したように、信頼度は、ｘnm,sの絶対値である。絶対値が大きくなるほど、信頼度は高くなっている。
　図１０は、本実施形態の応用例を示す説明図である。一例として、推定結果において信頼度が低い領域を除外して、何らかの解析などに使用することができる。図１０において、黒領域が、除外される領域である。
　さらに、位相の推定値に基づいて変位解析等が実行された後、解析結果を信頼できないような領域を可視表示するような用途にも、本実施形態による信頼度を応用可能である。
実施形態４．
　図１１は、第４の実施形態の画像解析装置の構成例を示すブロック図である。図１１に示す画像解析装置４０は、図７に示された画像解析装置３０の構成要素に加えて、前フィルタ部１６０をさらに含む。画像解析装置４０において、画素選択部１１０、低次元圧縮部１２０、展開部１３０、および、位相平均化部１４２と活性化部１４３とを含む最適化部１４０は、第３の実施形態におけるそれらと同様に動作する。ただし、本実施形態では、位相平均化部１４２は、第３の実施形態における位相平均化部１４２とは異なる評価関数を使用する。
　前フィルタ部１６０は、例えば、下記の（８）式によって、行列Γを計算する。（８）式において、"Ｗfilter"は、あらかじめ決められている重みである。"ｎ∈Ω（nm）は、選択されたＭ個の画素のうちのいずれかであることを示す。（８）式は、分子にｙn,sとｙn,s'の共役複素数との積を含むので、近傍観測画素の位相差をＷfilterで平均化することを意味する。（８）式で示される行列（以下、添え字が略された行列Γと表現する。）は、観測画像における何らかのノイズの特徴を表している。そして、行列Γは、本実施形態における位相平均化部１４２が使用する評価関数に反映される。
Figure JPOXMLDOC01-appb-M000005
　位相平均化部１４２は、図８におけるステップＳ１３３と同様の処理を実行する。しかし、本実施形態では、位相平均化部１４２は、以下の（９）式を用いて、ｘtempを算出する。その他の処理は、図７に示された画像解析装置３０の処理と同じである。（９）式は、（３）式におけるｙn,sとｙn,s'の共役複素数との積が、行列Γで置き換えられた式に相当する。
Figure JPOXMLDOC01-appb-M000006
　本実施形態では、前フィルタ部１６０が、近傍観測画素の位相差の平均が反映された重み行列を決定するので、観測画像内の情報から最適な重みを決めることができる。そして、位相平均化部１４２が、近傍観測画素の位相差の平均が反映された行列Γが含まれる評価関数を用いてｘtempを算出するので、位相ノイズ低減の程度をさらに大きくすることができる。
実施形態５．
　図１２は、第５の実施形態の画像解析装置の構成例を示すブロック図である。図１１に示す画像解析装置５０は、図１１に示された画像解析装置４０の構成要素に加えて、重み算出部１７０をさらに含む。
　重み算出部１７０は、下記の（１０）式で示すように、重みWnmを算出する。すなわち、重み算出部１７０は、重みWnmを、行列Γに基づいて算出する。
Figure JPOXMLDOC01-appb-M000007
　画像解析装置５０における他の動作は、第４の実施形態における動作と同じである。
　図１３は、ＣＰＵを有するコンピュータの一例を示すブロック図である。コンピュータは、上記の各実施形態の画像解析装置に実装される。ＣＰＵ１０００は、記憶装置１００１に格納されたプログラム（ソフトウェア要素：コード）に従って処理を実行することによって、上記の実施形態における各機能を実現する。すなわち、図１，図４，図５，図７，図１１，図１２に示された画像解析装置１０，１０Ａ，２０，３０，４０，５０における、画素選択部１１０、低次元圧縮部１２０、展開部１３０、最適化部１４０、位相平均化部１４２、活性化部１４３、平滑化部１５０、前フィルタ部１６０、および重み算出部１７０の機能を実現する。なお、ＣＰＵ１０００に代えて、ＧＰＵ、または、ＣＰＵとＧＰＵとの組み合わせを使用することもできる。
　記憶装置１００１は、例えば、非一時的なコンピュータ可読媒体（non-transitory computer readable medium ）である。非一時的なコンピュータ可読媒体は、様々なタイプの実体のある記録媒体（tangible storage medium ）を含む。非一時的なコンピュータ可読媒体の具体例として、磁気記録媒体（例えば、ハードディスク）、光磁気記録媒体（例えば、光磁気ディスク）、ＣＤ－ＲＯＭ（Compact Disc-Read Only Memory ）、ＣＤ－Ｒ（Compact Disc-Recordable ）、ＣＤ－Ｒ／Ｗ（Compact Disc-ReWritable ）、半導体メモリ（例えば、マスクＲＯＭ、ＰＲＯＭ（Programmable ROM）、ＥＰＲＯＭ（Erasable PROM ）、フラッシュＲＯＭ）がある。また、記憶装置１００１として、データを書き換え可能な記憶媒体が用いられる場合には、記憶装置１００１は、ＳＡＲ画像格納部２００および重み行列格納部３００としても使用可能である。また、記憶装置１００１は、保存部１３１，１４１としても使用可能である。
　また、プログラムは、様々なタイプの一時的なコンピュータ可読媒体（transitory computer readable medium ）に格納されてもよい。一時的なコンピュータ可読媒体には、例えば、有線通信路または無線通信路を介して、すなわち、電気信号、光信号または電磁波を介して、プログラムが供給される。
　メモリ１００２は、例えばＲＡＭ（Random Access Memory）で実現され、ＣＰＵ１０００が処理を実行するときに一時的にデータを格納する記憶手段である。メモリ１００２に、記憶装置１００１または一時的なコンピュータ可読媒体が保持するプログラムが転送され、ＣＰＵ１０００がメモリ１００２内のプログラムに基づいて処理を実行するような形態も想定しうる。
　図１４は、画像解析装置の主要部を示すブロック図である。図１４に示す画像解析装置１００は、同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する画素選択部（画素選択手段）１１（実施形態では、画素選択部１１０で実現される。）と、評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する低次元圧縮部（低次元圧縮手段）１２（実施形態では、低次元圧縮部１２０で実現される。）と、低次元圧縮部１２による圧縮結果を元の画素空間に戻し、空間相関位相推定値を算出する展開部（展開手段）１３（実施形態では、展開部１３０で実現される。）と、評価値を、空間相関位相推定値と、画素選択部１１によって選択された位置の画素とに近づけることによって、評価関数を最適化する最適化部（最適化手段）１４（実施形態では、最適化部１４０で実現される。）とを備えている。
　図１５は、他の態様の画像解析装置の主要部を示すブロック図である。図１５に示す画像解析装置１０１において、最適化部１４が使用する評価関数は、画素選択部１１によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして複数画像に亘って加算される項を含み、最適化部１４は、当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する制限部（制限手段）１５（実施形態では、活性化部１４３で実現される。）を含む。
　上記の実施形態の一部または全部は、以下の付記のようにも記載され得るが、以下に限定されるわけではない。
（付記１）同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する画素選択手段と、
　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する低次元圧縮手段と、
　前記低次元圧縮手段による圧縮結果を元の画素空間に戻し、空間相関位相推定値を算出する展開手段と、
　前記評価値を、前記空間相関位相推定値と、前記画素選択手段によって選択された位置の画素とに近づけることによって、前記評価関数を最適化する最適化手段と
　を備える画像解析装置。
（付記２）前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算される項を含み、
　前記最適化手段は、当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する制限手段を含む
　付記１の画像解析装置。
（付記３）前記低次元圧縮手段は、前記複素ベクトルに低次元圧縮するための行列を乗算し、
　前記展開手段は、前記低次元圧縮手段による圧縮結果に展開のための行列を乗算して前記空間相関位相推定値を導出する
　付記１または付記２の画像解析装置。
（付記４）前記低次元圧縮手段による圧縮結果を過去の圧縮結果で平滑化してから前記展開手段に供給する平滑化手段（実施形態では、平滑化部１５０で実現される。）をさらに備える
　付記１から付記３のうちのいずれかの画像解析装置。
（付記５）前記展開手段の出力の絶対値を保存する絶対値保存手段（実施形態では、保存部１３１で実現される。）をさらに備える
　付記１から付記４のうちのいずれかの画像解析装置。
（付記６）前記展開手段の出力の位相を保存する位相保存手段（実施形態では、保存部１３１で実現される。）をさらに備える
　付記１から付記５のうちのいずれかの画像解析装置。
（付記７）前記最適化手段の出力の絶対値を保存する絶対値保存手段（実施形態では、保存部１４１で実現される。）をさらに備える
　付記１から付記４のうちのいずれかの画像解析装置。
（付記８）前記最適化手段の出力の位相を保存する位相保存手段（実施形態では、保存部１４１で実現される。）をさらに備える
　付記１から付記４のうちのいずれかまたは付記７の画像解析装置。
（付記９）同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択し、
　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮し、
　前記低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出し、
　前記評価値を、前記空間相関位相推定値と、前記選択された位置の画素とに近づけることによって、前記評価関数を最適化する
　画像解析方法。
（付記１０）前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算された項を含み、
　当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する
　付記９の画像解析方法。
（付記１１）コンピュータに、
　同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する処理と、
　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する処理と、
　前記低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出する処理と、
　前記評価値を、前記空間相関位相推定値と、前記選択された位置の画素とに近づけることによって、前記評価関数を最適化する処理と
　を実行させる画像解析プログラムが格納されたコンピュータ読み取り可能な記録媒体。
（付記１２）前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算された項を含み、
　コンピュータに、当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する処理
　を実行させる画像解析プログラムが格納された付記１１のコンピュータ読み取り可能な記録媒体。
（付記１３）コンピュータに、
　同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する処理と、
　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する処理と、
　前記低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出する処理と、
　前記評価値を、前記空間相関位相推定値と、前記選択された位置の画素とに近づけることによって、前記評価関数を最適化する処理と
　を実行させるための画像解析プログラム。
（付記１４）前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算された項を含み、
　コンピュータに、当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する処理
　を実行させる付記１３の画像解析プログラム。
　以上、実施形態を参照して本願発明を説明したが、本願発明は上記の実施形態に限定されるものではない。本願発明の構成や詳細には、本願発明のスコープ内で当業者が理解し得る様々な変更をすることができる。
　１１　　　画素選択部
　１２　　　低次元圧縮部
　１３　　　展開部
　１４　　　最適化部
　１５　　　制限部
　１０，１０Ａ，２０，３０，４０，５０　画像解析装置
　１００，１０１　画像解析装置
　１１０　　画素選択部
　１２０　　低次元圧縮部
　１３１，１４１　保存部
　１３０　　展開部
　１４０　　最適化部
　１４２　　位相平均化部
　１４３　　活性化部
　１５０　　平滑化部
　１６０　　前フィルタ部
　１７０　　重み算出部
　２００　　ＳＡＲ画像格納部
　３００　　重み行列格納部
　１０００　ＣＰＵ
　１００１　記憶装置
　１００２　メモリ

Claims (12)

  1. 　同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する画素選択手段と、
    　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する低次元圧縮手段と、
    　前記低次元圧縮手段による圧縮結果を元の画素空間に戻し、空間相関位相推定値を算出する展開手段と、
    　前記評価値を、前記空間相関位相推定値と、前記画素選択手段によって選択された位置の画素とに近づけることによって、前記評価関数を最適化する最適化手段と
    　を備える画像解析装置。
  2. 　前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算される項を含み、
    　前記最適化手段は、当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する制限手段を含む
    　請求項１に記載の画像解析装置。
  3. 　前記低次元圧縮手段は、前記複素ベクトルに低次元圧縮するための行列を乗算し、
    　前記展開手段は、前記低次元圧縮手段による圧縮結果に展開のための行列を乗算して前記空間相関位相推定値を導出する
    　請求項１または請求項２に記載の画像解析装置。
  4. 　前記低次元圧縮手段による圧縮結果を過去の圧縮結果で平滑化してから前記展開手段に供給する平滑化手段をさらに備える
    　請求項１から請求項３のうちのいずれか１項に記載の画像解析装置。
  5. 　前記展開手段の出力の絶対値を保存する絶対値保存手段をさらに備える
    　請求項１から請求項４のうちのいずれか１項に記載の画像解析装置。
  6. 　前記展開手段の出力の位相を保存する位相保存手段をさらに備える
    　請求項１から請求項５のうちのいずれか１項に記載の画像解析装置。
  7. 　前記最適化手段の出力の絶対値を保存する絶対値保存手段をさらに備える
    　請求項１から請求項４のうちのいずれか１項に記載の画像解析装置。
  8. 　前記最適化手段の出力の位相を保存する位相保存手段をさらに備える
    　請求項１から請求項４のうちのいずれか１項または請求項７に記載の画像解析装置。
  9. 　同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択し、
    　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮し、
    　前記低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出し、
    　前記評価値を、前記空間相関位相推定値と、前記選択された位置の画素とに近づけることによって、前記評価関数を最適化する
    　画像解析方法。
  10. 　前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算された項を含み、
    　当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する
    　請求項９に記載の画像解析方法。
  11. 　コンピュータに、
    　同一地域が記録された複数画像のうちの１枚の画像における複数の位置の画素を選択する処理と、
    　評価関数を最適化したときの評価値としての複素ベクトルを低次元空間に低次元圧縮する処理と、
    　前記低次元圧縮された結果を元の画素空間に戻し、空間相関位相推定値を算出する処理と、
    　前記評価値を、前記空間相関位相推定値と、前記選択された位置の画素とに近づけることによって、前記評価関数を最適化する処理と
    　を実行させる画像解析プログラムが格納されたコンピュータ読み取り可能な記録媒体。
  12. 　前記評価関数は、前記画素選択手段によって選択された位置の画素の値を表す複素数とその共役複素数との積を重み付けして前記複数画像に亘って加算された項を含み、
    　コンピュータに、当該評価関数を最適化して得られる値を非線形変換して所定値以内に制限する処理
    　を実行させる画像解析プログラムが格納された請求項１１に記載のコンピュータ読み取り可能な記録媒体。
PCT/JP2020/036843 2020-09-29 2020-09-29 画像解析装置および画像解析方法 WO2022070242A1 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2022553249A JP7416277B2 (ja) 2020-09-29 2020-09-29 画像解析装置および画像解析方法
PCT/JP2020/036843 WO2022070242A1 (ja) 2020-09-29 2020-09-29 画像解析装置および画像解析方法
US18/027,687 US20230377201A1 (en) 2020-09-29 2020-09-29 Image analysis device and image analysis method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2020/036843 WO2022070242A1 (ja) 2020-09-29 2020-09-29 画像解析装置および画像解析方法

Publications (1)

Publication Number Publication Date
WO2022070242A1 true WO2022070242A1 (ja) 2022-04-07

Family

ID=80951490

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2020/036843 WO2022070242A1 (ja) 2020-09-29 2020-09-29 画像解析装置および画像解析方法

Country Status (3)

Country Link
US (1) US20230377201A1 (ja)
JP (1) JP7416277B2 (ja)
WO (1) WO2022070242A1 (ja)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012533051A (ja) * 2009-07-08 2012-12-20 テレ−リレヴァメント エウローパ−ティ．エルレ．エ． エッセ．エルレ．エルレ． 同一領域上で取得されたｓａｒ画像から得られるインターフェログラムのフィルタリング処理方法
CN103439708A (zh) * 2013-08-29 2013-12-11 西安电子科技大学 基于广义散射矢量的极化InSAR干涉图估计方法
JP2014095585A (ja) * 2012-11-08 2014-05-22 Mitsubishi Electric Corp 画像レーダ処理装置及び画像レーダ処理方法
WO2019106850A1 (ja) * 2017-12-01 2019-06-06 日本電気株式会社 Ｓａｒ画像解析システム、画像処理装置、画像処理方法および画像処理プログラム
US20190353779A1 (en) * 2018-05-15 2019-11-21 University Of Electronic Science And Technology Of China Ground-based interferometric synthetic aperture radar-based atmospheric phase compensation method

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2012533051A (ja) * 2009-07-08 2012-12-20 テレ−リレヴァメント エウローパ−ティ．エルレ．エ． エッセ．エルレ．エルレ． 同一領域上で取得されたｓａｒ画像から得られるインターフェログラムのフィルタリング処理方法
JP2014095585A (ja) * 2012-11-08 2014-05-22 Mitsubishi Electric Corp 画像レーダ処理装置及び画像レーダ処理方法
CN103439708A (zh) * 2013-08-29 2013-12-11 西安电子科技大学 基于广义散射矢量的极化InSAR干涉图估计方法
WO2019106850A1 (ja) * 2017-12-01 2019-06-06 日本電気株式会社 Ｓａｒ画像解析システム、画像処理装置、画像処理方法および画像処理プログラム
US20190353779A1 (en) * 2018-05-15 2019-11-21 University Of Electronic Science And Technology Of China Ground-based interferometric synthetic aperture radar-based atmospheric phase compensation method

Also Published As

Publication number Publication date
JPWO2022070242A1 (ja) 2022-04-07
JP7416277B2 (ja) 2024-01-17
US20230377201A1 (en) 2023-11-23

Similar Documents

Publication Publication Date Title
JP7315748B2 (ja) データ識別器訓練方法、データ識別器訓練装置、プログラム及び訓練方法
RU2503938C2 (ru) Способ и устройство частотного анализа данных
Yu et al. Forecasting of hydrologic time series with ridge regression in feature space
US11468273B2 (en) Systems and methods for detecting and classifying anomalous features in one-dimensional data
JP6727642B2 (ja) 学習型アルゴリズムによるフォーカス補正処理方法
US10339442B2 (en) Corrected mean-covariance RBMs and general high-order semi-RBMs for large-scale collaborative filtering and prediction
JP7020494B2 (ja) Ｓａｒ画像解析システム、画像処理装置、画像処理方法および画像処理プログラム
US20200293857A1 (en) Cnn processing device, cnn processing method, and program
US7526414B2 (en) Density morphing and mode propagation for Bayesian filtering
Brajard et al. Combining data assimilation and machine learning to emulate a dynamical model from sparse and noisy observations: a case study with the Lorenz 96 model
CN110275147B (zh) 基于迁移深度神经网络的人体行为微多普勒的分类识别方法
CN106483561B (zh) 一种频率域的子波分解方法
CN117111000A (zh) 一种基于双通道注意力残差网络的sar梳状谱干扰抑制方法
WO2022070242A1 (ja) 画像解析装置および画像解析方法
Ailliot et al. Dynamical partitioning of directional ocean wave spectra
CN113449637A (zh) 毫米波雷达人体骨骼姿势估计方法和装置
WO2022070243A1 (ja) 画像解析装置および画像解析方法
CN111856568A (zh) 基于mwv模型的频率域多道反射系数联合反演方法及系统
WO2021192038A1 (ja) 画像解析装置および画像解析方法
CN110223238B (zh) 编码光照成像重构方法及装置
WO2021186557A1 (ja) データ処理装置およびデータ処理方法
Abdechiri et al. Exploring chaotic attractors in nonlinear dynamical system under fractal theory
Hanias et al. On efficient multistep non-linear time series prediction in chaotic diode resonator circuits by optimizing the combination of non-linear time series analysis and neural networks
WO2022269754A1 (ja) 画像解析装置および画像解析方法
Mei-na et al. A Bayesian approach for jointly estimating the model order and the DOAs

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 20956172

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2022553249

Country of ref document: JP

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20956172

Country of ref document: EP

Kind code of ref document: A1