Elevation analysis device and elevation analysis method

Abstract

The present invention comprises: a DEM difference processing unit 14 that calculates a difference between interference SAR data and elevation data (DEM data) according to a numerical elevation model, to thereby acquire differential interference SAR data in which a displacement component and uncorrelated noise due to ground surface fluctuation are reduced; a filter processing unit 15 that performs high-pass filter processing on the differential interference SAR data, to thereby acquire phase difference data in which an atmospheric error and an orbital fringe error of a satellite are reduced; an elevation change calculation unit 16 that calculates elevation change data of the ground surface from the phase difference data acquired in this way; and an elevation calculation unit 18 that calculates elevation data by adding the elevation change data to the DEM data oversampled by an oversampling unit 17. The present invention makes it possible to accurately measure an elevation difference using interference SAR data according to a repeat pass observation technique, and accurately calculate the elevation of the ground surface from the elevation difference and the DEM data.

Description

Elevation analysis device and elevation analysis method

[0001] The present invention relates to an elevation analysis device and an elevation analysis method, and is particularly suitable for use in a technique for measuring the elevation of the ground surface using an interferometric SAR.

[0002] Elevation monitoring is a process of monitoring changes in terrain using a geographic information system (GIS) and remote sensing technology, and is utilized in a wide range of fields as follows. Elevation monitoring helps to take appropriate measures by grasping terrain changes in real time or collecting data periodically. ○ River basin flood control (grasping sediment movement and deposition volume, monitoring the status of sediment control dams and levees, etc.) ○ Measurement of sediment movement volume and range during sediment disasters ○ Monitoring of landfill ○ Measurement of terrain deformation during large-scale collapses such as dams, slopes, and mines

[0003] As techniques for performing elevation monitoring, there are ground surveying, airborne laser scanning, differential measurement using a digital elevation model (DEM), and variation measurement by interferometric SAR (InSAR) using satellite SAR (Synthetic Aperture Radar). Ground surveying and airborne laser scanning are relatively highly accurate, but they have the drawback that wide-area surveying is difficult and costly, so the survey frequency is limited. Measurement by DEM difference is low-cost when using an existing DEM that has already been measured and published, and can cover a wide area, but it is difficult to obtain timely data and is not suitable for real-time measurement.

[0004] While interferometric SAR (SAR) measurement of deformation can cover a wide area because it uses satellites and has high accuracy in detecting deformations in millimeters or centimeters, it cannot be used for phenomena such as those described above where the ground surface changes significantly. In contrast, a technique is known that uses SAR's temporal displacement (time-series differential interferometry) to calculate elevation (see, for example, Non-Patent Literature 1). Non-Patent Literature 1 discloses a method in which data imaged by SAR is analyzed to calculate the elevation of the target area in a planar manner, then image data after sediment movement is analyzed to calculate the planar height, and the amount of sediment movement is calculated by performing a differential analysis on these two planar height data.

[0005] However, performing time-series differential interferometry requires a large amount of SAR image data, which presents the problem of taking a long time to collect all the necessary data. Furthermore, since there should be no ground surface deformation due to phenomena other than those mentioned above within the time period of the time-series analysis, it is difficult to accurately calculate the amount of sediment transport in areas where deformation is likely to occur.

[0006] Furthermore, techniques have been proposed to detect surface deformation by combining SAR with DEM or Digital Terrain Model (DTM) (see, for example, Patent Documents 1 and 2). Patent Document 1 discloses the formation of N-1 differential interference images in relation to the main image for N (N>20) images acquired using SAR that can be used over many years, with a vertical accuracy better than 50 meters. The differential interference images are created by subtracting the influence of the terrain from the phase difference of each pixel and using an existing DEM.

[0007] Furthermore, Patent Document 1 explains that when a satellite is equipped with two radar antennas, one for transmitting and one for receiving, that illuminate the same surface area, and they are located at distances ρ and ρ+Δρ from a point on the ground, the measured value of the phase difference Φq with respect to the elevation q of that point on the ground can be extrapolated with an accuracy of several meters using the following equation (a). Here, λ is the wavelength of the radar wave, B is the baseline length which is the distance between two satellites projected perpendicular to the line of sight, and θ is the angle of incidence of the radar wave emitted from the satellite, which is shifted from the direction perpendicular to the ground surface.

[0008]

[0009] Furthermore, Patent Document 1 explains that, after creating a temporary series of interferometer phases for each pixel of a differential interference image selected based on the statistical properties of reflectance, the elevation difference Δq at a single point can be calculated by the following equation (b) by calculating the phase linear component a (the slope a of the straight line in the equation Φ = aB + C (where C is a constant) for the minimum square) for each temporary series with respect to the baseline.

[0010]

[0011] In the terrain deformation extraction device described in Patent Document 2, interference fringe images are generated by processing two sets of SAR reconstructed images acquired from one artificial satellite. At the same time, DTM is generated from the interference fringe images and two sets of orbital data, and a pseudo-interference fringe image is generated by using the DTM and two sets of orbital data for the same region, assuming that there is no terrain deformation. Then, the amount of terrain deformation is calculated by processing the interference fringe image and pseudo-interference fringe image generated in this way, and then geometric distortion is calculated from the amount of terrain deformation and DTM and a geometric transformation is performed.

[0012] Special Publication No. 2003-500658 Japanese Patent Application Publication No. 8-15426

[0013] November 1, 2022, "NEC and others develop technology to monitor elevation from SAR data - for managing sediment transport in rivers," Internet search [URL: https: / / uchubiz.com / article / new9328 / ]

[0014] There are two techniques for measuring elevation using interferometric SAR: the bistatic observation technique, which measures the difference between SAR images taken simultaneously by two satellites at the same location, and the repeat-pass observation technique, which measures the difference between SAR images taken twice by a single satellite at the same location. The former allows for high-precision measurements but requires a dedicated satellite. On the other hand, the latter does not require a dedicated satellite, but because the two observations are not simultaneous, atmospheric errors, ground deformation, and reduced coherence significantly reduce the measurement accuracy.

[0015] Patent Document 1 discloses that the elevation difference Δq can be calculated as shown in equation (b), but it does not disclose how to reduce the effects of atmospheric errors, ground deformation, and reduced coherence when using repeat-pass observation technology. Therefore, it is difficult to accurately measure the elevation difference using the technology described in Patent Document 1.

[0016] Furthermore, Patent Document 1 discloses a technique for measuring motion in urban areas and landslide zones, but does not disclose how to calculate elevation. On the other hand, Patent Document 2 discloses generating a DTM from interference fringe images and satellite orbit data, but does not disclose a method for accurately calculating elevation values ​​while reducing the effects of atmospheric errors, ground deformation, and reduced coherence.

[0017] This invention was made to solve such problems, and aims to enable accurate measurement of the elevation at the time of satellite observation in areas where elevation differences may occur due to phenomena in which the ground surface changes significantly in a relatively short period of time, using interferometric SAR data obtained by repeat pass observation technology.

[0018] To solve the above-mentioned problems, the present invention calculates the difference between interference SAR data generated by a predetermined interference process on data observed by a synthetic aperture radar mounted on a satellite and elevation data from a numerical elevation model. This calculates differential interference SAR data from which displacement components due to ground deformation and uncorrelated noise have been reduced. Then, by applying a high-pass filter to the differential interference SAR data, phase difference data for each unit area is obtained from the differential interference SAR data from which atmospheric errors and satellite orbital fringe errors have been further reduced. From the phase difference data thus obtained, elevation change data for each unit area of ​​the ground is calculated using known information regarding ground surface observation by synthetic aperture radar. Furthermore, elevation data from a numerical elevation model is oversampled to a spatial resolution corresponding to the size of the unit area of ​​elevation change data, and elevation change data for each unit area is calculated by adding the elevation change data to the oversampled elevation data for each unit area.

[0019] While it is possible to calculate topography from interferometric SAR data generated based on observation data from synthetic aperture radar, when interferometric SAR data is generated using repeat-pass observation technology, various noises such as atmospheric errors, ground deformation, satellite orbital fringe errors, and uncorrelated noise are included in the interferometric SAR data due to changes in conditions during two observations of the same location. In contrast, according to the present invention configured as described above, phase difference data with the above-mentioned noises reduced is generated from the interferometric SAR data, and elevation change data is calculated from this phase difference data using known information. Therefore, using interferometric SAR data obtained by repeat-pass observation technology, elevation differences related to phenomena in which the ground surface changes significantly in a relatively short period of time can be measured with high accuracy. Furthermore, by adding the elevation change data calculated with high accuracy in this way to elevation data that has been oversampled to match the spatial resolution of the elevation change data, elevation data for each unit area of ​​the ground surface can be calculated with high accuracy. As a result, according to the present invention, it is possible to accurately measure the elevation at the time of satellite observation in areas where elevation differences may occur due to phenomena in which the ground surface changes significantly in a relatively short period of time, using interferometric SAR data obtained by repeat pass observation technology.

[0020] This figure shows an example of the network configuration of the analysis system according to this embodiment. This is a block diagram showing an example of the functional configuration of the elevation analysis device according to the first embodiment. This figure schematically shows the relationship between elevation change and baseline length that can be detected within the range where the dynamic range of the phase difference is 2π or less. This figure shows a state in which the sensitivity of CSK is partially expanded. This figure schematically shows the baseline length for multiple SAR images obtained when multiple observations are made using TSX. This is a block diagram showing an example of the functional configuration of the elevation analysis device according to the second embodiment. This figure shows an example of data processed by the mask processing unit and the elevation change calculation unit. This figure shows an example of data processed by the interpolation processing unit. This is a block diagram showing an example of the functional configuration of the elevation analysis device according to the third embodiment. This figure is for explaining layover and shadow. This is a block diagram showing an example of the functional configuration of the elevation analysis device according to the fourth embodiment.

[0021] Hereinafter, one embodiment of the present invention will be described with reference to the drawings. Figure 1 is a diagram showing an example of the network configuration of the analysis system 100 according to this embodiment. The analysis system 100 of this embodiment includes a satellite base station 50 that receives satellite data from one satellite S flying in a satellite orbit in space, a satellite database 51 that records satellite data, a DEM database 52 that stores elevation data (for example, DEM data) based on a digital elevation model, and an elevation analysis device 10 that analyzes the satellite data recorded in the satellite database 51.

[0022] Satellite S is equipped with a synthetic aperture radar (SAR) and photographs a point on Earth from a predetermined orbit, and transmits satellite data, including the captured ground image (hereinafter referred to as SAR image or SAR data), to the satellite base station 50. The SAR photographs the same point multiple times at different times, and each time transmits satellite data, including the SAR image, to the satellite base station 50.

[0023] The satellite base station 50 receives satellite data from satellite S and records it in the satellite database 51. The satellite database 51 is connected to a communication network N such as the Internet. The DEM database 52 is also connected to the communication network N. The DEM database 52 stores pre-generated DEM data.

[0024] The elevation analysis device 10 acquires satellite data stored in the satellite database 51 and DEM data pre-stored in the DEM database 52 via the communication network N. Using this acquired data, it performs processing related to calculating elevation differences (elevation changes) caused by phenomena that cause large changes in the ground surface over a relatively short period of time, and the elevation at the time of satellite observation. In this embodiment, the phenomena targeted for elevation difference calculation are, for example, watershed flood control, sediment disasters, embankments, and large-scale collapses, which are phenomena in which elevation changes of approximately 1 m to several tens of m occur over a short period of time.

[0025] Figure 2 is a block diagram showing an example of the functional configuration of the elevation analysis device 10 according to the first embodiment. As shown in Figure 2, the elevation analysis device 10 according to the first embodiment includes, as a functional configuration, a satellite data acquisition unit 11, an interferometric SAR image generation unit 12, a DEM data acquisition unit 13, a DEM difference processing unit 14, a filter processing unit 15, an elevation change calculation unit 16, an oversampling unit 17, and an elevation calculation unit 18.

[0026] The above-described functional blocks 11 to 18 perform the processes described below through the cooperation of hardware and software. For example, the processes of the above-described functional blocks 11 to 18 are executed by the operation of programs stored in storage media such as RAM, ROM, hard disk, or semiconductor memory, under the control of a microcomputer configured to include a CPU, RAM, ROM, etc. In addition to the microcomputer, a DSP (Digital Signal Processor) may also be included.

[0027] The physical configuration described herein is illustrative and does not necessarily have to be independent. For example, the elevation analysis device 10 may include an LSI (Large-Scale Integration) that integrates the CPU with RAM and ROM. Furthermore, although this embodiment describes the case where the elevation analysis device 10 is composed of a single computer, the elevation analysis device 10 may be realized by combining multiple computers.

[0028] The satellite data acquisition unit 11 acquires satellite data stored in the satellite database 51. In this embodiment, repeat pass observation technology is used, and the satellite data acquisition unit 11 acquires satellite data from two sets of images from among multiple sets of satellite data generated by multiple SAR imagings for the region of interest for which elevation difference is calculated.

[0029] The satellite data acquired by the satellite data acquisition unit 11 includes, in addition to SAR images, metadata containing known information regarding surface observations by SAR. This metadata includes, for example, radar wave wavelength λ and baseline length B. ⊥This includes the off-nadia angle θ, inclination distance r, satellite position, satellite velocity, satellite direction of movement (satellite azimuth angle), satellite time, and geographical information of the observation area. Baseline length B ⊥ is the distance between observation points in two observations by one satellite S. The off-nadia angle θ is the angle of incidence of the radar waves emitted from the SAR, deviating from the direction perpendicular to the Earth's surface. The slope distance r is the distance from the antenna of satellite S to the observation point on the Earth's surface.

[0030] The interferometric SAR image generation unit 12 generates an interferometric SAR image (interferometric SAR data) by performing a predetermined interference process to interfere two SAR images acquired by the satellite data acquisition unit 11. Known processes related to interferometric SAR can be used as the predetermined interference process.

[0031] Interferometric SAR (SAR) is a technique that uses the phase of radar waves to obtain images. It involves conducting two observations of the same location on the Earth's surface to generate two SAR images (phase images), and then interfering them and taking the difference to obtain information about even small distance differences as phase difference. The image generated by interfering these two SAR images is called an interferometric SAR image.

[0032] While typical interferometric SAR can capture surface movement occurring during two observation periods as a change in distance between the satellite and the surface in millimeters or centimeters, this embodiment is not intended to capture such minute surface displacements. The interferometric SAR image generation unit 12 simply generates an interferometric SAR image by interfering two SAR images, and the phase difference (φ, described later) corresponding to the minute surface displacements mentioned above is sufficient. disp It is not necessary to calculate the elevation difference. However, the generated interferometric SAR image can be said to be a phase difference image that includes phase difference information associated with ground deformation, which is different from the phenomenon targeted for elevation difference calculation.

[0033] The following equation (1) shows multiple elements included in the phase difference component Δφ of the interference SAR image generated by the interference SAR image generation unit 12. atm This is atmospheric error, φ orb φ is the orbital fringe error of satellite S. disp φ is the displacement component due to ground surface deformation.topo is the DEM differential phase, φ noise is the uncorrelated noise, and 2kπ is the phase integer value bias.

[0034]

[0035] When the atmosphere changes during two observations of the SAR image, the delay state of the radar wave phase changes due to the changed atmospheric state. Therefore, even if there is no change on the ground surface at the two observation time points, if there is a change in the atmospheric state, a phase difference will occur between the two SAR images accordingly. This phase difference is included in the phase difference component Δφ of the interferometric SAR image as the atmospheric error φ atm and comes to be included in the phase difference component Δφ of the interferometric SAR image.

[0036] Also, since the orbit of satellite S does not become exactly the same during two observations of the SAR image, even if there is no change on the ground surface at the two observation time points, a phase difference will occur between the two SAR images due to the difference in orbits. This appears as striped noise in the interferometric SAR image. Since the orbit of satellite S is known, the orbit fringe error φ orb can be estimated based on the information of its orbit and removed to some extent. However, it is difficult to completely remove it, and the orbit fringe error φ orb comes to be included in the phase difference component Δφ of the interferometric SAR image.

[0037] The displacement component φ due to ground surface movement disp is the information we want to obtain in the case of a general interferometric SAR as described above, but in the case of this embodiment, this is treated as noise information. On the other hand, the information we want to obtain in this embodiment is the DEM differential phase φ topo which is. The DEM differential phase φ topo is information indicating the phase difference from the terrain shown by DEM data as will be described later. The uncorrelated noise φ noise is noise that occurs randomly.

[0038] The phase integer bias of 2kπ can be described as an error related to the so-called phase unwrapping problem. In interferometric SAR, the phase value obtained is the fractional part when the slope distance r from the satellite S antenna to the observation point on the Earth's surface is divided by the wavelength λ of the radar wave. Therefore, the obtained phase value is folded (wrapped) into the range of {-π, π}. For this reason, for phases with a large dynamic range, there is uncertainty in integer multiples of 2π. That is, it is not known how many cycles k of the radar wave 2π occurred, only that it was at a specific position within the last cycle. It is possible to estimate the value of cycle k by a predetermined unwrapping process, but it is not always possible to estimate it accurately, and this is included as noise in the phase difference component Δφ of the interferometric SAR image.

[0039] As shown in equation (1), the phase difference component Δφ of the interference SAR image generated by the interference SAR image generation unit 12 includes the desired DEM difference phase φ. topo In addition, atmospheric error φ atm , orbital fringe error φ orb , displacement component φ due to ground surface deformation disp , uncorrelated noise φ noise This also includes a phase integer bias of 2kπ. In this embodiment, as detailed below, by performing processing to remove or reduce these noises, the change in elevation of the ground surface can be accurately calculated from the phase difference component Δφ' from which various noises have been reduced, and furthermore, the elevation value of the ground surface can be accurately calculated from this change in elevation.

[0040] The DEM data acquisition unit 13 acquires DEM data that is pre-stored in the DEM database 52. The DEM data acquired here is elevation data of the region of interest for calculating the elevation difference, and includes data of the same region as the region of interest included in the SAR data acquired by the satellite data acquisition unit 11. It is not necessary for the region of the SAR data acquired by the satellite data acquisition unit 11 and the region of the DEM data acquired by the DEM data acquisition unit 13 to be exactly the same; it is sufficient that the region of interest is included in both sets of data. Geographic information is attached as metadata to both the SAR data and the DEM data, and for example, it is possible to identify the region of interest using this geographic information.

[0041] The DEM difference processing unit 14 corresponds to the difference processing unit in the claims, and calculates the difference between the interference SAR data generated by the interference SAR image generation unit 12 and the DEM data acquired by the DEM data acquisition unit 13, thereby calculating the displacement component φ due to ground deformation. disp and uncorrelated noise φ noise The DEM difference processing unit 14 obtains differential interferometric SAR data (hereinafter sometimes referred to as differential interferometric SAR image) with reduced interference. At this time, the DEM difference processing unit 14 geocodes the interferometric SAR data and the DEM data to align them and then calculates the difference.

[0042] DEM data is existing elevation data of a numerical model generated by aerial laser surveying or other methods at a specific point in the past. For example, it is data that represents the terrain by arranging elevation values ​​in a grid. In contrast, interferometric SAR data is data that represents the slight displacement of the ground surface at two observation points in two SAR images acquired from the satellite database 51 as a phase difference for each unit area (each pixel). Since the phase difference is information corresponding to the distance difference, it is possible to obtain differential interferometric SAR data corresponding to the distance difference from the existing elevation data by calculating the difference between the interferometric SAR data and the DEM data.

[0043] Here, by shortening the time interval between the two observations in which two SAR images are taken, it is possible to reduce the possibility of ground deformation occurring during that time. In this case, the differential interferometric SAR data generated by the DEM differential processing unit 14 is the displacement component φ due to ground deformation. disp It can be said that the displacement component φ is reduced. This differential interference SAR data is disp In addition, uncorrelated noise φ noise While this noise has been reduced, other noises remain.

[0044] The filter processing unit 15 performs a high-pass filter process on the differential interference SAR data (differential interference SAR image) generated by the DEM difference processing unit 14, thereby further reducing atmospheric error φ from the differential interference SAR data. atm and the orbital fringe error φ of satellite S orb Phase difference data for each unit region with reduced distortion is obtained. In particular, the filter processing unit 15 of this embodiment performs low-pass filtering on the differential interference SAR data and subtracts the data obtained therefrom from the original differential interference SAR data before the low-pass filtering.

[0045] In other words, the filter processing unit 15 first extracts smooth change components from the differential interferometric SAR image by performing a low-pass filter on the differential interferometric SAR image. Then, it removes these extracted smooth change components from the original differential interferometric SAR image. By applying a high-pass filter through this process, atmospheric errors φ, which have components in the low-frequency region, are removed. atm and orbital fringe error φ orb It is possible to obtain phase difference data (differential interferometric SAR images) with significantly reduced distortion.

[0046] Note that atmospheric error φ atm orbital fringe error φ orb Since these phenomena occur on the order of several hundred meters to several kilometers, it is preferable to perform filtering with a window size that matches their size. With a window size of this size, it will be larger than the region of interest for which the elevation difference is calculated, and therefore will not significantly affect the phase information of the region of interest to be calculated.

[0047] The elevation change calculation unit 16 calculates elevation change data for each unit area of ​​the Earth's surface from the phase difference data generated by the filter processing unit 15, using known information regarding Earth's surface observation by SAR. The known information is information included as metadata in the satellite data acquired by the satellite data acquisition unit 11, and includes the wavelength λ of the radar wave and the baseline length B. ⊥ The off-nadia angle θ and the slope distance r are used. The elevation change calculation unit 16 uses this known information to calculate the elevation change dh for each unit area using the following equation (2).

[0048]

[0049] As a result of the processing up to the filter processing unit 15, in the above equation (1), the atmospheric error φ included in the phase difference component Δφ of the interferometric SAR image is atm , orbital fringe error φ orb , displacement component φ due to ground surface deformation disp and uncorrelated noise φ noise This can be significantly reduced. As for the remaining integer phase bias 2kπ, by performing observations by SAR under conditions such that the phase difference data generated by the filter processing unit 15 is within the range where phase unwrapping processing is unnecessary, the difficult-to-estimate phase unwrapping processing is avoided, and the integer phase bias 2kπ can be ignored. This makes it possible to use equation (2) above.

[0050] For this purpose, in this embodiment, the baseline length B is as shown in formula (2). ⊥ The sensitivity of elevation change dh, which depends on the phase difference, was theoretically calculated for each satellite S, and observation conditions were sought under which phase unwrapping is unnecessary when detecting elevation change dh of magnitude related to the phenomenon of the elevation difference calculation target. Here, the sensitivity of elevation change dh means that elevation change dh of magnitude related to the phenomenon of the elevation difference calculation target can be detected within the range where the dynamic range of the phase difference is 2π or less (not exceeding {-π, π}).

[0051] Figure 3 shows the elevation change dh and baseline length B that can be detected within the range where the dynamic range of the phase difference is 2π or less. ⊥This diagram schematically shows the relationship. In Figure 3, the horizontal axis is the baseline length B. ⊥ The vertical axis represents the elevation change dh, and the range in which elevation change dh can be detected within the range where the dynamic range of the phase difference is 2π or less is shown for each of several types of satellites. Range 301 shows the sensitivity when satellite S is Sentinel-1 in the C band (wavelength λ = 5.6 cm), range 302 shows the sensitivity when satellite S is ALOS-2 in the L band (wavelength λ = 24 cm), range 303 shows the sensitivity when satellite S is TSX in the X band (wavelength λ = 3.1 cm), and range 304 shows the sensitivity when satellite S is CSK in the X band (wavelength λ = 3.1 cm).

[0052] In Figure 3, the minimum value of the elevation change dh that can be detected within the range where the dynamic range of the phase difference is 2π or less is dh. min and maximum value dh max This was calculated as shown in equation (3), which was derived from equation (2). The minimum value of the elevation change dh shown in equation (3) min This is based on the case of low coherence. The minimum value dh when coherence is high. min This is calculated using equation (4).

[0053]

[0054] Based on the sensitivity calculation results shown in Figure 3, satellite S and baseline length B are selected to have sensitivity to the amount of elevation change related to the elevation difference calculation target phenomenon within the range where the dynamic range of the phase difference is 2π or less. ⊥ The system is configured to perform SAR observations under the condition that the combination of the above conditions is met. The DEM difference processing unit 14 then generates differential interferometric SAR data using the interferometric SAR data generated from the data observed under those conditions.

[0055] In Figure 3, the sensitivity of Sentinel-1, shown in range 301, is given by baseline length B ⊥ This indicates that in the range of 100m or less, elevation changes dh of approximately 12m or less are undetectable (baseline length B ⊥ The sensitivity at 100m is dh min = about 12m, dh max(= approximately 57 m). For example, if the amount of elevation change related to the phenomenon of the elevation difference calculation target is only a few meters, Sentinel-1 is unsuitable (does not meet the conditions) because it does not have sensitivity to such elevation changes dh of a few meters. Sentinel-1 has a usable baseline length B ⊥ From the perspective of having a narrow scope, it cannot be said to be desirable.

[0056] The sensitivity of ALOS-2 shown in range 302 is given by baseline length B ⊥ This indicates that elevation changes dh of approximately 10 m or less are undetectable within a range of 500 m or less (baseline length B ⊥ The sensitivity at 500m is dh min = about 10m, dh max (= approximately 44 m). For example, if the amount of elevation change related to the phenomenon of the elevation difference calculation target is only a few meters, ALOS-2 would also be unsuitable because it does not have sensitivity to such elevation changes dh of a few meters.

[0057] The sensitivity of the TSX shown in range 303 is given by baseline length B ⊥ This indicates that elevation changes dh of approximately 2m or less are undetectable within a range of 300m or less (baseline length B ⊥ The sensitivity at 300m is dh min = approximately 2m, dh max (= approximately 8 m). For example, if the amount of elevation change related to the phenomenon of the elevation difference calculation target is around a few meters (but 2 meters or more), then TSX is appropriate (meets the conditions) because it is sensitive to such elevation changes of a few meters dh.

[0058] The sensitivity of the CSK shown in range 304 is given by baseline length B ⊥ This indicates that in the range between 600m and a predetermined value (1000m, not shown in the figure), elevation changes of approximately 2m or less are undetectable (baseline length B). ⊥ The sensitivity at 600m is dh min = approximately 2m, dh max (= approximately 5 m). For example, if the amount of elevation change related to the phenomenon of the elevation difference calculation target is around a few meters (but 2 meters or more), then CSK would also be appropriate because it is sensitive to such elevation changes dh of a few meters.

[0059] From the above, it can be seen that when the amount of elevation change related to the phenomenon of the elevation difference calculation target is on the order of a few meters, it is appropriate to perform observations with the TSX or CSK X-band satellite S. That is, using the X-band satellite S, the baseline length B shown in Figure 3 is appropriate. ⊥ If SAR observations are performed within this range, it is possible to detect the elevation change dh of the elevation difference calculation target within a range where the dynamic range of the phase difference is 2π or less, and the phase integer value bias of 2kπ can be ignored.

[0060] Figure 4 shows a magnified view of the CSK sensitivity shown in Figure 3. When the elevation change related to the elevation difference calculation target is 2m or more and 8m or less, the baseline length B ⊥ If the base length is less than 100m, the sensitivity of the elevation change dh is too high, while the baseline length B ⊥ When the baseline length is 300m or more, the sensitivity of the elevation change dh tends to be insufficient. Therefore, when using the above elevation change amount as the target for elevation difference calculation using CSK, the baseline length B ⊥ It is preferable to conduct SAR observations in a range of 100m to 300m. The same principle applies when using TSX.

[0061] Furthermore, if the amount of elevation change of the elevation difference calculation target changes, the appropriate satellite S and baseline length B will be selected accordingly. ⊥ Needless to say, the combination will change.

[0062] One of the things necessary to satisfy the condition that the system is sensitive to the amount of elevation change of the elevation difference calculation target within a range where the dynamic range of the phase difference is 2π or less is an appropriate baseline length B. ⊥ One approach is to use two SAR images. In other words, it is important to obtain appropriate satellite data for two observations from the satellite data for multiple observations stored in the satellite database 51.

[0063] Figure 5 shows the baseline length B for multiple SAR images obtained when multiple observations were performed using TSX. ⊥This diagram schematically illustrates the situation. In Figure 5, the horizontal axis represents the time axis, the vertical axis represents the distance from the reference orbit to the actual observation position, and each circle indicates the distance at each observation point. The difference in the distance between two circles represents the baseline length B for the two observations. ⊥ This indicates.

[0064] In this embodiment, the satellite data acquisition unit 11 selects from among multiple satellite data corresponding to the multiple circles shown in Figure 5, the baseline length B ⊥ This selects and acquires a pair of satellite data that satisfies the above conditions. Using this as the first condition, in this embodiment, as described above, another pair of satellite data is acquired that further satisfies the second condition regarding a short time interval between the two observations. When the time interval between the two observations is short, the possibility of ground deformation occurring during that time is low, so the displacement component φ of ground deformation is added to the phase difference component Δφ of the interferometric SAR image. disp It is possible to exclude it or minimize it.

[0065] For example, if there are multiple combinations of satellite data that satisfy the first condition, one pair of satellite data with the shortest time interval may be selected from among them. Alternatively, if there are multiple combinations of satellite data that satisfy the first condition, one pair of satellite data whose time interval is less than the threshold may be selected from among them.

[0066] The oversampling unit 17 oversamples the DEM data acquired from the DEM database 52 by the DEM data acquisition unit 13 to a spatial resolution corresponding to the size of the unit area (pixel) of the elevation change data calculated by the elevation change calculation unit 16.

[0067] As described above, DEM data is data that represents terrain by arranging elevation values ​​in a grid, and generally available DEM data has a spatial resolution of 30m. In contrast, the spatial resolution of the elevation change data calculated by the elevation change calculation unit 16 is only a few meters. For example, the oversampling unit 17 oversamples the DEM data so that the spatial resolution of the DEM data matches the spatial resolution of the elevation change data.

[0068] The elevation calculation unit 18 calculates elevation data for each unit area by adding the elevation change data calculated by the elevation change calculation unit 16 to the DEM data that has been oversampled by the oversampling unit 17, for each unit area. The elevation change data calculated by the elevation change calculation unit 16 is data that shows the difference from the DEM data at the time of observation of the satellite data acquired by the satellite data acquisition unit 11. Therefore, the elevation data obtained by adding this elevation change data to the DEM data represents elevation values ​​that reflect the topography at the time of observation.

[0069] As explained in detail above, in the first embodiment, the difference between interferometric SAR data and elevation data from a digital elevation model is calculated to determine the displacement component φ due to ground deformation. disp and uncorrelated noise φ noise After obtaining differential interference SAR data with reduced interference, a high-pass filter is applied to the differential interference SAR data to further reduce atmospheric error φ from the differential interference SAR data. atm and the orbital fringe error φ of satellite S orb The phase difference data Δφ' (=φ) has been reduced. topo The phase difference data Δφ' obtained in this way is used to calculate the elevation change dh of the ground surface using equation (2) with known information regarding ground surface observation by SAR. In addition, the DEM data from the numerical elevation model is oversampled to a spatial resolution corresponding to the size of the unit area of ​​the elevation change data, and the elevation change dh data is added to the oversampled DEM data to calculate the elevation data.

[0070] According to the first embodiment configured in this way, phase difference data Δφ' with the various noises mentioned above reduced is generated from the interferometric SAR data, and elevation change dh is calculated from this phase difference data Δφ' using known information. Therefore, using interferometric SAR data obtained by repeat pass observation technology, elevation change dh related to phenomena in which the ground surface changes significantly, for example in units of meters, over a relatively short period of time can be measured with high accuracy. Furthermore, by adding the elevation change data calculated in this high accuracy to DEM data that has been oversampled to match the spatial resolution of the elevation change data, the elevation data of the ground surface can be calculated with high accuracy.

[0071] Furthermore, in the first embodiment, the baseline length B is sensitive to the amount of elevation change of the elevation difference calculation target within a range where the dynamic range of the phase difference is 2π or less. ⊥ By using satellite data observed under these conditions to acquire differential interferometric SAR data, it is possible to eliminate the need for phase unwrapping, which is generally associated with estimation errors, and to effectively ignore the phase integer value bias of 2kπ. This improves the detection accuracy of elevation change dh, and thus improves the detection accuracy of the elevation value calculated based on said elevation change dh.

[0072] Figure 6 is a block diagram showing an example of the functional configuration of the elevation analysis device 10A according to the second embodiment. In Figure 6, components with the same reference numerals as those shown in Figure 2 have the same function, so redundant explanations are omitted here.

[0073] As shown in Figure 6, the elevation analysis device 10A according to the second embodiment has an elevation change calculation unit 16A and an elevation calculation unit 18A instead of the elevation change calculation unit 16 and elevation calculation unit 18 in Figure 2, and further includes a mask processing unit 21 and an interpolation processing unit 22. The processing of these functional blocks 16A, 18A, 21, and 22 is also executed by the operation of a program stored in a storage medium under the control of a microcomputer.

[0074] The mask processing unit 21 removes phase difference data from the phase difference data for each unit area filtered by the filter processing unit 15, specifically from areas where the spatial dispersion of the phase difference is greater than a threshold. The elevation change calculation unit 16A calculates elevation change data for areas where the phase difference data has not been removed by the mask processing unit 21.

[0075] For example, the mask processing unit 21, using the phase difference data for each unit region filtered by the filter processing unit 15, moves a predetermined-sized window while shifting it, calculates the variance of the phase difference group contained within each window region, and removes the phase difference data for window regions where the variance is greater than a threshold. As an example, the variance calculated for each window region is normalized to a value between 0 and 1, and the phase difference data for window regions where the normalized variance is greater than 0.5 is removed. The threshold may be made adjustable by the user.

[0076] Figure 7 shows an example of data processed by the mask processing unit 21 and the elevation change calculation unit 16A. Figure 7(a) shows the phase difference data (phase difference image representing the phase difference for each pixel) filtered by the filter processing unit 15. This Figure 7(a) shows the phase difference data for a part of the region of interest for which elevation is calculated, and there is a variation region 701 with a large spatial dispersion of phase differences around the middle of this region.

[0077] Figure 7(b) shows the phase difference data obtained by processing the phase difference data in Figure 7(a) with the mask processing unit 21. The area 702, which appears to be filled in black, is the area where the phase difference data has been removed because the spatial dispersion of the phase difference is greater than the threshold. Figure 7(c) shows the elevation change data (elevation change image showing elevation change dh for each pixel) calculated from the masked phase difference data as in Figure 7(b). Elevation change data has not been calculated for the masked variation area 703.

[0078] As described above, by performing the processing of the mask processing unit 21, it is possible to calculate elevation change data by excluding low-reliability unit-domain phase difference data that has the properties of uncorrelated noise. For example, it is possible to remove phase difference data in areas where the spatial dispersion of phase difference is large, such as the water surface, and calculate elevation change data for the other areas.

[0079] The interpolation processing unit 22 uses the elevation change data calculated by the elevation change calculation unit 16A for the areas where data was not removed by the mask processing unit 21, and generates elevation change data corresponding to the areas where data was removed by the mask processing unit 21 through interpolation.

[0080] As shown in Figure 7(c), the processing performed by the mask processing unit 21 and the elevation change calculation unit 16A results in missing data in the generated elevation change data. Before the elevation calculation unit 18A adds the elevation change data to the DEM data, the interpolation processing unit 22 generates the elevation change data for the missing regions by interpolation. The interpolation process can be performed by applying, for example, a Gaussian filter, Laplace interpolation, nearest neighbor interpolation, bilinear interpolation, bicubic interpolation, spline interpolation, Lanczos interpolation, or other known image interpolation filters.

[0081] The elevation calculation unit 18A calculates elevation data by adding elevation change data calculated by the elevation change calculation unit 16A for areas where data was not removed by the mask processing unit 21 and elevation change data generated by the interpolation processing unit 22 for areas where data was removed by the mask processing unit 21 to the DEM that has been oversampled by the oversampling unit 17.

[0082] Figure 8 shows an example of data processed by the interpolation processing unit 22. Figure 8(a) shows elevation change data generated by the elevation change calculation unit 16A from phase difference data masked by the mask processing unit 21. Figure 8(a) shows elevation change data for a portion of the region of interest for which elevation is calculated, and the area that appears to be blacked out within the dotted line-enclosed area 801 indicates a missing area where there is no elevation change data.

[0083] Figure 8(b) shows elevation change data obtained by performing interpolation processing on the missing region in Figure 8(a) using the interpolation processing unit 22. The elevation calculation unit 18A calculates elevation data by adding the elevation change data generated as shown in Figure 8(b) to the DEM data oversampled by the oversampling unit 17. In this way, elevation data without missing regions can be obtained.

[0084] Furthermore, elevation data calculated in areas where elevation change data has been interpolated is less reliable than elevation data calculated in areas where it has not been interpolated. Therefore, it is necessary to recognize the interpolated areas when performing geometric correction or interpretation of elevation data, for example. To this end, the elevation calculation unit 18A may generate information that can identify the interpolated areas and record it in association with the elevation data.

[0085] Figure 9 is a block diagram showing an example of the functional configuration of the elevation analysis device 10B according to the third embodiment. In Figure 9, components with the same reference numerals as those shown in Figure 6 have the same function, so redundant explanations are omitted here.

[0086] As shown in Figure 9, the elevation analysis device 10B according to the third embodiment further includes an invisible area estimation unit 31 in addition to the configuration shown in Figure 6. The processing of this invisible area estimation unit 31 is also executed by the operation of a program stored in a storage medium, controlled by a microcomputer.

[0087] In this example, we show a configuration that further includes an invisible area estimation unit 31 in addition to the functional blocks of the second embodiment shown in Figure 6, but a configuration that further includes an invisible area estimation unit 31 in addition to the functional blocks of the first embodiment shown in Figure 2 is also possible.

[0088] The invisible area estimation unit 31 estimates the area of ​​elevation data corresponding to the invisible area with respect to at least one of layover and shadow, based on the satellite orbit data acquired by the satellite data acquisition unit 11 and the elevation data calculated for each unit area by the elevation calculation unit 18A.

[0089] Figure 10 is a diagram illustrating layover and shadow. Figure 10(a) is an explanatory diagram of layover. As shown in Figure 10(a), on slopes where the terrain gradient is steeper than the incidence angle α (> off-nadia angle θ) of the radar waves emitted from satellite S to the ground surface, the distance to satellite S is shorter at a higher elevation location B (ridge) located further back than at a lower elevation location A (valley) located closer to the foreground. Since satellite S records satellite data in the order in which it receives reflected waves from the ground surface, the higher elevation location B is recorded before the lower elevation location A, resulting in a phenomenon where the reception order is reversed. As a result, the interferometric SAR image is inverted vertically, and the phase difference on the near side of location B is not correctly observed and disappears. This is the invisible region due to layover. Layover is more likely to occur when the incidence angle α is small.

[0090] Figure 10(b) is an explanatory diagram of shadows. As shown in Figure 10(b), on slopes where the terrain gradient is steeper than the supplemental angle β of the incident angle of radar waves at the Earth's surface, the radar waves on the opposite slope (the back side of the ridge) as seen from satellite S are blocked by the ridge and cannot be observed, resulting in a loss of data in this area. This is the invisible region due to shadows. Shadows tend to occur more easily when the supplemental angle β of the incident angle is small.

[0091] Invisible areas due to layover and invisible areas due to shadow can be estimated by simulation using information such as the off-nadia angle θ of radar waves emitted from satellite S, the azimuth angle of satellite S, and topographic information of the Earth's surface. There are known techniques for estimating invisible areas from DEM data.

[0092] The invisible area estimation unit 31 estimates the invisible area from the entire range of elevation data by performing a simulation using known technology on the elevation data calculated by the elevation calculation unit 18A, based on the metadata included in the satellite data acquired by the satellite data acquisition unit 11 and the elevation data calculated by the elevation calculation unit 18A.

[0093] The invisible regions in the elevation data calculated by the elevation calculation unit 18A are different from the invisible regions in the DEM data acquired by the DEM data acquisition unit 13. Therefore, in order to recognize these invisible regions when interpreting the elevation data, for example, the invisible region estimation unit 31 generates information that can identify the estimated invisible regions and records it in association with the invisible regions of the elevation data.

[0094] Figure 11 is a block diagram showing an example of the functional configuration of the elevation analysis device 10C according to the fourth embodiment. In Figure 11, components with the same reference numerals as those shown in Figure 9 have the same function, so redundant explanations are omitted here.

[0095] As shown in Figure 11, the elevation analysis device 10C according to the fourth embodiment further includes a second filter processing unit 41 in addition to the configuration shown in Figure 9. The processing of this second filter processing unit 41 is also executed by the operation of a program stored in a storage medium under the control of a microcomputer.

[0096] In this example, we show a configuration that further includes a second filter processing unit 41 in addition to the functional block of the third embodiment shown in Figure 9. However, the configuration may also include a second filter processing unit 41 in addition to the functional block of the first embodiment shown in Figure 2, or a second filter processing unit 41 in addition to the functional block of the second embodiment shown in Figure 6.

[0097] The second filter processing unit 41 obtains phase difference data with reduced uncorrelated noise by performing a filter process on the phase difference data for each unit region filtered by the filter processing unit 15 to remove specific frequency components. For example, an adaptive filter processing unit is provided as the second filter processing unit 41 to further perform a process to remove specific frequency components. This removes uncorrelated noise φ that may remain even after processing by the interference SAR image generation unit 12. noiseThis can be further reduced. For example, if a SAR image includes areas of water such as the sea or a lake, interference in the SAR image is less likely to occur in those areas, resulting in a grainy noise. This can be reduced by adaptive filtering.

[0098] The mask processing unit 21 removes phase difference data from the phase difference data filtered by the second filter processing unit 41 in regions where the spatial variance of the phase difference is greater than a threshold.

[0099] Although the first to fourth embodiments have been described above, the present invention is not limited thereto. In some cases, SAR observations may be performed under conditions including a phase integer bias of 2kπ, and the elevation change dh may be calculated from equation (2) after performing phase unwrapping. For example, a satellite S and baseline length B that completely cover the amount of elevation change targeted for elevation difference calculation within a range of 2π ⊥ If no such combination exists, satellite S and baseline length B cover the elevation change amount targeted for elevation difference calculation within a range of 4π. ⊥ Alternatively, SAR observations can be performed in combination with the above, and the elevation change dh can be calculated from equation (2) after performing phase unwrapping. In this case as well, since the range estimated by phase unwrapping is narrower, the occurrence of estimation errors can be suppressed compared to when normal phase unwrapping is performed.

[0100] Furthermore, in the first to fourth embodiments described above, a Gaussian filter may be applied to the image representing the elevation change to remove fine noise after the elevation change calculation units 16 and 16A.

[0101] Furthermore, the theoretical uncertainty (error) of the elevation change data calculated by the elevation change calculation units 16 and 16A and the elevation data calculated by the elevation calculation units 18 and 18A may also be calculated. Equation (5) below is an error index value σ that indicates the theoretical uncertainty of the elevation change data. dh This is the formula for calculating [the value].

[0102]

[0103] In Equation (5), σ γ represents the standard deviation of the phase associated with the decrease in coherence γ, and is calculated for each pixel from the value of coherence γ. σ atm represents the standard deviation of the phase due to the residual atmospheric error after the high-pass filter processing by the filter processing unit 15, and one value is calculated for the entire filtered phase difference image.

[0104] The error index value σ dh is, σ γ and σ atm and represents the standard deviation of the elevation difference considering the sensitivity of the elevation change dh, and is calculated for each pixel. In Equation (5), the part of dh 2π / 2π corresponds to the operation for converting the standard deviation of the phase into the standard deviation of the elevation difference.

[0105] The following Equation (6) is an equation for calculating the error index value σ DEMref indicating the theoretical uncertainty of the elevation data. In Equation (6), σ DEMex represents the standard deviation of the DEM data acquired by the DEM data acquisition unit 13, and is calculated for each pixel. The error index value σ DEMref represents the standard deviation of the elevation data, and is calculated for each pixel.

[0106]

[0107] Examples of configurations applicable to the elevation analysis device according to the embodiment are summarized below.

[0108] [Configuration 1] An elevation analysis device comprising: a difference processing unit that calculates the difference between interference SAR data generated by a predetermined interference processing on data observed by a synthetic aperture radar mounted on a satellite and elevation data from a numerical elevation model, thereby acquiring differential interference SAR data in which displacement components due to ground deformation and uncorrelated noise have been reduced; a filter processing unit that performs a high-pass filter processing on the differential interference SAR data to acquire phase difference data for each unit area from the differential interference SAR data in which atmospheric errors and orbital fringe errors of the satellite have been further reduced; an elevation change calculation unit that calculates elevation change data for each unit area of ​​the ground surface from the phase difference data using known information regarding ground surface observation by the synthetic aperture radar; an oversample unit that oversamples the elevation data from the numerical elevation model to a spatial resolution corresponding to the size of the unit area of ​​the elevation change data; and an elevation calculation unit that calculates elevation data for each unit area by adding the elevation change data calculated by the elevation change calculation unit to the oversampled elevation data for each unit area.

[0109] [Configuration 2] The elevation analysis device according to Configuration 1, further comprising a mask processing unit that removes phase difference data for regions where the spatial dispersion of the phase difference is greater than a threshold in the phase difference data for each unit region filtered by the filter processing unit, wherein the elevation change calculation unit calculates the elevation change data for regions where the phase difference data has not been removed by the mask processing unit.

[0110] [Configuration 3] An elevation analysis device according to Configuration 2, further comprising an interpolation processing unit that generates elevation change data corresponding to the areas from which data was removed by the mask processing unit by interpolation using the elevation change data calculated by the elevation change calculation unit for the areas from which data was not removed by the mask processing unit, wherein the elevation calculation unit calculates the elevation data by adding the elevation change data calculated by the elevation change calculation unit for the areas from which data was not removed by the mask processing unit and the elevation change data generated by the interpolation processing unit for the areas from which data was removed by the mask processing unit to the oversampled elevation data.

[0111] [Configuration 4] An elevation analysis device according to any one of Configurations 1 to 3, further comprising an invisible area estimation unit that estimates an area of ​​elevation data corresponding to an invisible area with respect to at least one of layover and shadow, based on the orbital data of the satellite and the elevation data calculated for each unit area by the elevation calculation unit.

[0112] [Configuration 5] An elevation analysis device according to any one of Configurations 1 to 4, further comprising a second filter processing unit that performs a filter process to remove specific frequency components from the phase difference data for each unit region filtered by the above filter processing unit, thereby obtaining phase difference data with reduced uncorrelated noise.

[0113] [Configuration 6] An elevation analysis device according to any one of Configurations 2 to 4, further comprising a second filter processing unit that performs a filter process to remove specific frequency components from the phase difference data for each unit region filtered by the above filter processing unit, thereby obtaining phase difference data with reduced uncorrelated noise, wherein the mask processing unit removes phase difference data in the phase difference data filtered by the second filter processing unit in regions where the spatial dispersion of the phase difference is greater than a threshold. [Configuration 7] An elevation analysis device according to any one of Configurations 1 to 6, wherein the known information includes the baseline length, which is the distance between observation points in two observations by one satellite, and the difference processing unit obtains the difference interferometric SAR data using data observed under conditions where the baseline length is sensitive to the amount of elevation change of the elevation difference calculation target within a range where the dynamic range of the phase difference shown in the phase difference data is 2π or less. [Configuration 8] The elevation analysis device according to any one of Configurations 1 to 6, wherein the above known information includes the baseline length, which is the distance between observation points in two observations by one satellite, and the difference processing unit acquires the difference interferometric SAR data using two observation data that satisfy a first condition in which the baseline length is sensitive to the amount of elevation change of the elevation difference calculation target within a range in which the dynamic range of the phase difference shown in the above phase difference data is 2π or less, and a second condition in which the time interval between the two observations is short.

[0114] Furthermore, the first to fourth embodiments described above are merely examples of how the present invention may be implemented, and the technical scope of the present invention should not be interpreted as being limited by them. In other words, the present invention can be implemented in various forms without departing from its gist or its main features.

[0115] 10, 10A, 10B, 10C Elevation analysis device 11 Satellite data acquisition unit 12 Interferometric SAR image generation unit 13 DEM data acquisition unit 14 DEM difference processing unit (difference processing unit) 15 Filter processing unit 16, 16A Elevation change calculation unit 17 Oversampling unit 18, 18A Elevation calculation unit 21 Mask processing unit 22 Interpolation processing unit 31 Invisible area estimation unit 41 Second filter processing unit 51 Satellite database 52 DEM database S Satellite

Claims

1. An elevation analysis device comprising: a difference processing unit that calculates the difference between interference SAR data generated by a predetermined interference processing on data observed by a synthetic aperture radar mounted on a satellite and elevation data from a numerical elevation model, thereby acquiring differential interference SAR data with reduced displacement components due to ground deformation and uncorrelated noise; a filter processing unit that performs a high-pass filter processing on the differential interference SAR data to acquire phase difference data for each unit area from the differential interference SAR data, thereby further reducing atmospheric errors and orbital fringe errors of the satellite; an elevation change calculation unit that uses known information regarding ground surface observation by the synthetic aperture radar to calculate elevation change data for each unit area of ​​the ground surface from the phase difference data; an oversample unit that oversamples the elevation data from the numerical elevation model to a spatial resolution corresponding to the size of the unit area of ​​the elevation change data; and an elevation calculation unit that calculates elevation data for each unit area by adding the elevation change data calculated by the elevation change calculation unit to the oversampled elevation data for each unit area.

2. The elevation analysis device according to claim 1, further comprising a masking unit that removes phase difference data for regions where the spatial dispersion of the phase difference is greater than a threshold in the phase difference data for each unit region filtered by the above-mentioned filtering unit, wherein the elevation change calculation unit calculates the elevation change data for regions where the phase difference data has not been removed by the above-mentioned masking unit.

3. The elevation analysis device according to claim 2, further comprising an interpolation processing unit that generates elevation change data corresponding to the areas from which data was removed by the mask processing unit by interpolation using the elevation change data calculated by the elevation change calculation unit for the areas from which data was not removed by the mask processing unit, wherein the elevation calculation unit calculates the elevation data by adding the elevation change data calculated by the elevation change calculation unit for the areas from which data was not removed by the mask processing unit and the elevation change data generated by the interpolation processing unit for the areas from which data was removed by the mask processing unit to the oversampled elevation data.

4. An elevation analysis device according to any one of claims 1 to 3, further comprising an invisible area estimation unit that estimates an area of ​​elevation data corresponding to an invisible area with respect to at least one of layover and shadow, based on the orbital data of the satellite and the elevation data calculated for each unit area by the elevation calculation unit.

5. The elevation analysis apparatus according to claim 2 or 3, further comprising a second filter processing unit that performs a filter process to remove specific frequency components from the phase difference data for each unit region filtered by the above filter processing unit, thereby obtaining phase difference data with reduced uncorrelated noise, wherein the mask processing unit removes phase difference data in regions where the spatial variance of the phase difference is greater than a threshold in the phase difference data filtered by the second filter processing unit.

6. The elevation analysis device according to claim 1, wherein the above known information includes the baseline length, which is the distance between observation points in two observations by one satellite, and the difference processing unit acquires the difference interferometric SAR data using data observed under conditions of a baseline length that is sensitive to the amount of elevation change of the elevation difference calculation target within a range in which the dynamic range of the phase difference shown in the above phase difference data is 2π or less.

7. The elevation analysis device according to claim 1, wherein the above known information includes the baseline length, which is the distance between observation points in two observations by one satellite, and the difference processing unit acquires the difference interferometric SAR data using two observation data that satisfy a first condition, which is that the baseline length is sensitive to the amount of elevation change of the elevation difference calculation target within a range where the dynamic range of the phase difference shown in the phase difference data is 2π or less, and a second condition, which is that the time interval between the two observations is short.

8. The computer's differential processing unit calculates the difference between interference SAR data generated by a predetermined interference process on data observed by a synthetic aperture radar mounted on a satellite and elevation data from a digital elevation model, thereby obtaining differential interference SAR data with reduced displacement components due to ground deformation and uncorrelated noise; the computer's filter processing unit performs a high-pass filter process on the differential interference SAR data to obtain phase difference data for each unit area from the differential interference SAR data, thereby further reducing atmospheric errors and satellite orbital fringe errors; the computer's elevation change calculation unit calculates elevation change data for each unit area of ​​the ground surface from the phase difference data using known information regarding ground surface observation by the synthetic aperture radar; and the computer's oversampling unit oversamples the elevation data from the digital elevation model to a spatial resolution corresponding to the size of the unit area of ​​the elevation change data. An elevation analysis method characterized by comprising the step of the elevation calculation unit of the above-mentioned computer calculating elevation data for each unit area by adding the elevation change data calculated by the elevation change calculation unit to the oversampled elevation data for each unit area.

Images