Glacier ice reserve estimation method based on flight type ice measuring radar

Abstract

The invention relates to a glacier ice reserve estimation method based on a flight type ice measuring radar. The method comprises the following steps: acquiring a large-range terrain of a glacier; emitting microwaves for measuring the ice thickness in the glacier boundary; receiving echoes for measuring the ice thickness in the glacier boundary; positioning the acquired microwave signals in real time; processing ice thickness data; and estimating the glacier reserves. In order to improve the glacier detection accuracy, the core equipment glacier thickness measuring radar has centimeter-level ice layer distinguishing capability and five hundred-meter-level ice thickness penetrating capability, has the characteristics of light weight, miniaturization, high reliability and the like, is suitable for high-cold and high-altitude conditions, and is suitable for being carried by an unmanned aerial vehicle platform. Through the ultra-wideband high-stability low-phase noise signal generation technology, the weak echo signal receiving technology, the ultra-wideband antenna technology and the like, the thickness measuring radar system has the effects of light weight, miniaturization, high reliability and the like. The unmanned aerial vehicle is suitable for climatic conditions such as high glacier altitude and thin air, the defect that manual arrival or pull-type measurement must be relied on is overcome, and the testing efficiency is improved.

Description

Technical Field

[0001] The invention relates to a glacier ice reserve estimation method based on a flying ice measuring radar, which is a new hydrological method for calculating glacier water storage. Background Art

[0002] There are three main methods for estimating glacier reserves. The first is the empirical formula method, including the area-volume empirical formula. The second is the ice thickness model estimation method. A variety of physical models of varying complexity have emerged to simulate ice thickness distribution, such as those based on glacier material conservation and ice flow mechanics principles, shallow ice approximate theory, glacier flow laws, glacier velocity method or minimum method. These estimation models infer ice thickness from glacier surface characteristics, but whether the model parameters can be transplanted to glaciers in other regions has not been proven, and the uncertainty of the parameters is unknown, resulting in large uncertainties in the model estimation results. The third is direct field measurement methods, such as the early drilling method, flower pole and snow pit method, which are the most direct methods for measuring glacier thickness, but these methods are time-consuming and labor-intensive, and can only obtain very sparse point thickness information.

[0003] The three methods have their own advantages and disadvantages. In general, the first method is only suitable for the rough estimation of regional glacier reserves. When applied to a single glacier, the accuracy can only be kept at the same order of magnitude as the actual result, and there are large differences in the volume-area formula parameters of different types of glaciers. The second method is only suitable for the estimation of glacier reserves in a small area, and there is a lack of sufficient measured data for model parameter calibration. The third method of field measurement is almost impossible to achieve, because the plateau climate conditions where the glacier is located are harsh, such as high altitude and thin air, and precise measurement is difficult and inefficient. How to quickly, efficiently and accurately measure and estimate glacier water storage is a key issue that needs to be solved. Summary of the invention

[0004] In order to overcome the shortcomings of the prior art, the present invention proposes a method for estimating glacier ice reserves based on a flying ice radar. The method proposes a method for using a drone equipped with a glacier thickness radar to accurately measure the ice thickness of glaciers over a large area, accurately calculating the glacier water storage based on the measured glacier distribution and ice thickness, and providing an accurate scientific basis for the utilization and development of glacier water resources.

[0005] The object of the present invention is achieved as follows: a method for estimating glacier ice reserves based on a flying ice-measuring radar, wherein the system used in the method comprises: an unmanned aerial vehicle and a glacier thickness-measuring radar mounted on the unmanned aerial vehicle, and a control computer arranged on the ground; the glacier thickness-measuring radar comprises an ultra-wideband antenna unit, a transmitter, a receiver, an RTK positioning device capable of achieving centimeter-level positioning accuracy, and a signal processing device; the ultra-wideband antenna unit is a linear dipole transmitting antenna and a receiving antenna; the radar pulse signal of the transmitter is generated by using an avalanche triode and a step recovery diode to form a differential structure and a Marx cascade structure to form a carrier-free pulse generating circuit, and adopts a multi-tube joint working mode; the receiver is provided with a sensitivity time control circuit and a time-varying gain amplifier circuit to enhance the reflected signal of the target and increase the detection depth of the radar;

[0006] The steps of the method are as follows:

[0007] Step 1: Obtaining large-scale glacier topography: Using public satellite remote sensing data or using topographic mapping equipment for on-site mapping, obtain the boundary and surface elevation data of the glacier being studied;

[0008] Step 2, transmitting radar waves for measuring ice thickness within the boundary of the glacier: planning a survey line within the boundary of the glacier, and setting a flight path in the UAV according to the survey line; starting the UAV to cruise according to the pre-set flight path, turning on the transmitter of the glacier thickness radar, using avalanche triodes and SRD devices to generate carrier-free pulses that meet the requirements of radar detection of glaciers, and using a multi-tube cascade working mode to generate ultra-wideband low phase noise pulse transmission signals, and transmitting detection radar waves through a linear dipole antenna;

[0009] Step 3, receiving and measuring the radar echo reflected in the glacier: the receiver uses a sensitivity time control circuit and a time-varying gain amplifier circuit to increase the amplifier gain of the receiver in a nonlinear form as the pulse propagation time increases, so as to enhance the reflected signal of the target and increase the detection depth of the radar; the processing of the received signal adopts an equivalent sampling method to reconstruct the radar echo signal;

[0010] Step 4: Real-time positioning of the collected radar wave signal: In the process of receiving the echo, synchronize RTK to obtain the corresponding position coordinates. After the RTK information is collected, the longitude and latitude information can be used to track and dynamically display the measurement trajectory on the electronic map. At the same time, the longitude and latitude information is stored in the database table corresponding to the measurement data. Combined with GIS, the location information corresponding to the measurement data can be displayed;

[0011] Step 5, ice thickness data processing: The processing includes the following sub-steps:

[0012] Sub-step 5.1, filtering and gain adjustment, improving radar echo quality and measurement accuracy;

[0013] Sub-step 5.2, convert the two-way travel time of the radar signal into the ice thickness H to form the ice thickness layer on the survey line. The calculation formula is:

[0014]

[0015] Where: c is the propagation speed of electromagnetic waves in vacuum; t is the two-way travel time of electromagnetic waves in the measurement medium; ε is the relative dielectric constant of the medium; v is the propagation speed of radar pulses in the medium; d is the distance between the radar transmitting antenna and the receiving antenna;

[0016] Sub-step 5.3, the obtained ice thickness data is marked on the survey line, and the data includes ice thickness information and corresponding geographical location coordinates to form a point layer;

[0017] Sub-step 5.4, superimposing the ice thickness point layer to the glacier boundary, superimposing the thickness point layer data on the glacier map according to the geographic coordinates, so that the thickness information on the ice layer measurement line is added to the glacier boundary;

[0018] Sub-step 5.5, interpolate to obtain the ice thickness distribution map of all locations within the glacier boundary; in actual measurement, only the ice thickness at the survey line location can be obtained. Using the linear interpolation method, based on the existing ice thickness data on the survey line, the ice thickness of the area between the survey lines is estimated, thereby obtaining the ice thickness distribution data of all locations within the glacier boundary;

[0019] Step 6, glacier reserve estimation: Ice reserve estimation specifically includes the following sub-steps:

[0020] Sub-step 6.1, based on the ice thickness distribution data on the glacier radar survey line, form an ice thickness contour distribution map; divide it into 5m or 10m intervals along the elevation direction, and draw a map in GIS;

[0021] Sub-step 6.2, estimate the bottom area occupied by the glacier body at the same thickness in the GIS map;

[0022] Sub-step 6.3, the volume of the glacier block is obtained by multiplying the bottom area by the partition thickness; the volumes of the blocks at other ice thickness heights can be obtained in turn; and the total volume of the measured glacier body can be obtained by cumulative summation.

[0023] The advantages and beneficial effects of the present invention are as follows: the present invention uses a drone equipped with an ice-penetrating radar to directly measure the ice thickness distribution of a large range of glaciers by low-altitude flight, and then obtains the glacier ice reserves in combination with the glacier bedrock terrain, adapts to the requirements of harsh environmental conditions such as low temperature and vibration, and accurately calculates the distribution and thickness of the glacier, with the characteristics of high efficiency, high precision, high reliability, miniaturization, and wide adaptability. In order to improve the accuracy of glacier detection, the core equipment of the present invention, the glacier thickness radar, has an ice layer resolution capability of the order of centimeters and an ice thickness penetration capability of the order of 500 meters. The use of drones is suitable for climatic conditions such as high altitude and thin air of glaciers, overcoming the shortcomings of relying on manpower to reach or towed measurements in the past, reducing risks, and significantly improving the test efficiency, with the characteristics of fast, large-scale, and high-precision measurement. The ice thickness data measured by the present invention is accurate, and ice reserves can be efficiently calculated, with great engineering promotion value. BRIEF DESCRIPTION OF THE DRAWINGS

[0024] The present invention will be further described below in conjunction with the accompanying drawings and embodiments.

[0025] Figure 1 is a system principle block diagram used in the method described in the embodiment of the present invention;

[0026] Figure 2 is a flow chart of the method described in an embodiment of the present invention;

[0027] Figure 3 It is a schematic diagram of the basic principle of glacier thickness radar measurement according to an embodiment of the present invention;

[0028] Figure 4 It is a schematic diagram of equivalent sampling of the improved glacier thickness measuring radar described in an embodiment of the present invention. DETAILED DESCRIPTION

[0029] Example:

[0030] This embodiment is a method for estimating glacier ice reserves based on a flying ice radar. The system used in the method is as follows: Figure 1 As shown, it includes: a UAV 1, a glacier thickness measuring radar carried on the UAV, and a control computer 2 located on the ground. Figure 1The glacier thickness radar includes an ultra-wideband antenna unit, a transmitter, a receiver, an RTK (Real Time Kinematic) positioning device capable of achieving centimeter-level positioning accuracy, and a signal processing device; the ultra-wideband antenna unit is a linear dipole transmitting antenna and a receiving antenna; the radar pulse signal of the transmitter is generated by using avalanche triodes and step recovery diodes to form a differential structure and a Marx cascade structure to form a carrier-free pulse generation circuit, and a multi-tube joint working mode is adopted; the receiver is provided with a sensitivity time control circuit and a time-varying gain amplifier circuit, that is, the gain of the receiver amplifier increases in a nonlinear form (mainly an exponential form) as the pulse propagation time increases, which is used to enhance the target's reflected signal and increase the radar's detection depth. The receiver uses an equivalent sampling method to reconstruct the radar echo signal, and the receiver and the control computer are wirelessly connected.

[0031] The core equipment of the glacier thickness radar is carried by a drone. The drone can adopt a 4-axis 8-propeller multi-rotor configuration, with a maximum load of 30 kg and a maximum range of 16 kilometers when fully loaded. Conventional fixed-wing drones can also be used. The whole machine can adapt to the working environment temperature of -20℃ to 40℃, with a maximum flight altitude of 8,000 meters. It integrates large load, long range, strong signal, and high wind resistance, and is suitable for transportation in high-cold and high-altitude and various emergency scenarios.

[0032] The technical indicators of the receiver include: sensitivity: 1mV, sampling rate is 2GHz (sampling interval 0.5ns).

[0033] Ultra-wideband antenna unit: In order to avoid the influence of clutter, in areas with complex electromagnetic environments, glacier thickness radars generally need to use shielded butterfly antennas. The above antennas have the advantage of good radar echo signals, but they are large in size and weight. Since the glacier is in an uninhabited area with little clutter, in order to achieve lightweight and miniaturization, this embodiment adopts the following measures:

[0034] Based on theoretical simulation and experimental measurement, the glacier thickness radar uses a wideband radar wave with a center frequency of 100MHz. Considering the use of drones, in order to reduce the weight of the antenna, an impedance loading method based on Wu-King theory is adopted to broaden the bandwidth of the antenna through distributed resistance loading, so that it has better time domain and wideband characteristics. After the glacier thickness radar selects the above-mentioned unshielded antenna, the antenna arm length is 75cm, the total length is 150cm, the width is 10cm, and the weight is about 2kg (excluding the transmitter and receiver).

[0035] The steps of the method are as follows: Figure 2 As shown:

[0036] Step 1: Obtain large-scale topography of the glacier: Use public satellite remote sensing data or use topographic mapping equipment for on-site mapping to obtain the boundary of the glacier being studied and the glacier surface elevation data; obtaining these topographic data and surface elevation data is a prerequisite for calculating the thickness of the glacier.

[0037] The terrain and elevation data of the earth's surface can be obtained from the public data or models on the Internet. For the large-scale terrain of the glacier, the remote sensing image data (such as Landsat 8, source: USGS of the United States) and the digital surface model (such as DSM, source: JAXA of Japan) terrain data can be used to obtain the glacier boundary. The above-mentioned drone can also be used to directly carry out terrain mapping to obtain the glacier boundary and glacier surface elevation data.

[0038] Step 2, transmitting radar waves for measuring ice thickness within the boundary of the glacier: planning a survey line within the boundary of the glacier, and setting a flight path in the UAV according to the survey line; starting the UAV to cruise according to the pre-set flight path, turning on the transmitter of the glacier thickness radar, using avalanche triodes and SRD devices to generate carrier-free pulses that meet the requirements of radar detection of glaciers, and using a multi-tube cascade working mode to generate ultra-wideband low phase noise pulse transmission signals, and transmitting detection radar waves through a linear dipole antenna;

[0039] The generation of high-quality ultra-wideband pulse transmission signals is the core technology of the detection radar of this embodiment. The high voltage and wide bandwidth of the pulse signal can improve the detection distance and resolution of the radar, and the stability of its waveform has an impact on the imaging quality. The existing tunnel diodes and step recovery diodes cannot meet the requirements of the pulse source of the ground penetrating radar. Therefore, this embodiment uses an avalanche tube as a high-speed switch of the pulse source, which can generate very high narrow pulses.

[0040] Glacier measurement radar has strict requirements on the upper limit of repetition rate and stability of pulses of high-resolution detection radar pulse transmitters. These indicators directly affect the accuracy and detection distance of radar detection. Therefore, in order to meet the application requirements of pulse transmitters of high-resolution radars, the traditional MARX pulse generation circuit is debugged and improved, and the triggering form of the circuit is changed from series triggering to parallel triggering, so that the upper limit of repetition rate and pulse stability are greatly improved, the maximum repetition rate can reach 400KHz, and the waveform jitter is very small. The improved nanosecond ground penetrating radar transmitter circuit includes a power supply circuit, a trigger signal module, a drive amplifier module, and an avalanche transistor pulse generation circuit.

[0041] The basic principle of glacier thickness radar measurement is shown in Figure 3The transmitting antenna carried by the drone transmits electromagnetic waves in a direction toward the ice surface and penetrates the ice layer. The high-frequency electromagnetic waves are reflected by the sub-ice strata or target bodies with electrical differences and then return to be received by the receiving antenna. The two-way travel time t of the electromagnetic waves in the measurement medium can be calculated using the following formula:

[0042]

[0043] Where: H is the thickness of the measured medium; v is the propagation speed of the radar pulse in the medium; d is the distance between the radar transmitting antenna and the receiving antenna; the calculation formula of v is:

[0044]

[0045] Where: c is the propagation speed of electromagnetic waves in vacuum (0.29979m / ns), ε is the relative dielectric constant of the medium.

[0046] Therefore, the thickness of the measured medium can be calculated using the following formula:

[0047]

[0048] The difference in dielectric properties of different materials is a prerequisite for ground penetrating radar to detect target media. The dielectric constant of air is 1, ice is between 3-5, and rock (silt) is between 5-30. The dielectric properties of ice and the strata under the ice are significantly different, which provides a basis for accurate measurement of ground penetrating radar. When the high-frequency electromagnetic waves emitted by ground penetrating radar propagate, the electromagnetic field intensity and waveform will change with the dielectric properties and geometric shape of the medium. The radar waves received by the receiving antenna are analyzed for amplitude, waveform and frequency to obtain radar images, and the thickness of the target medium can be determined by image post-processing.

[0049] Step 3, receiving and measuring the radar echo reflected in the glacier: According to radar theory and practice, in the process of propagation of a carrier-free pulse signal penetrating ice and water, the signal energy decays according to an exponential law as the penetration depth increases. In order to better compensate for the nonlinear attenuation, the receiver uses a sensitivity time control (STC) circuit and a time-varying gain (TVG) amplifier circuit, so that the amplifier gain of the receiver increases in a nonlinear form (mainly an exponential form) as the pulse propagation time increases, which is used to enhance the reflected signal of the target and increase the detection depth of the radar; the processing of the received signal adopts an equivalent sampling method to reconstruct the radar echo signal;

[0050] Equivalent sampling is a transformation sampling that utilizes the periodic or quasi-periodic repetition characteristics of the received signal and performs sampling only once or several times in each period. After multiple repetitive periods of the signal, all samples that can reconstruct the waveform of the signal within one period are obtained.

[0051] The traditional sequential equivalent sampling scheme of ground penetrating radar only collects one point in each time period. If the information of 512 points of the echo signal is to be collected, 512 time periods will be required. For the high-speed flying ice-measuring radar, it is difficult to meet the system's requirements for acquisition speed and acquisition lateral resolution. In addition, if the antenna moves too fast, it is difficult to ensure the similarity of the 512 echo signals, which will result in the reconstructed signal being unable to be used for actual measurement. This embodiment improves the equivalent sampling to reduce the number of equivalent sampling periods. The principle is as follows:

[0052] Assume that the repetition period is T and the sampling time within the period is t p , collect K points, and the interval between sampling points in the cycle is τ Δ , the original signal can be reconstructed by collecting M cycles in total. The equivalent sampling step time between each cycle is Δt. A total of N points are collected, and the relationship between them satisfies:

[0053] K = T / τ Δ (4)

[0054] M = τ Δ / Δt (5)

[0055] N = M × K (6)

[0056] Sampling time t mk It can be expressed as:

[0057] t mk =(T+Δt)+ k τ Δ = mT s + k τ Δ (7)

[0058] Where: k = 0, 1, 2, ..., K-1; m = 0, 1, 2, ..., M-1. The traditional equivalent only collects one point each time, and this embodiment collects K points each time. From formula (5), it can be seen that the equivalent sampling cycle number M depends on the interval τ between the sampling points in the cycle Δ And the step time Δt between cycles, because Δt also determines the accuracy of equivalent sampling, so without sacrificing accuracy, Δt is reduced to reduce the number of cycles and improve the equivalent sampling speed. In this embodiment, Δt is the sampling interval of the sampling chip, so the system working rate depends on the working frequency of the AD chip. The improved equivalent sampling working process is as follows Figure 4 shown.

[0059] Figure 4 This is a schematic diagram of the hybrid equivalent sampling operation after the sequential equivalent sampling is improved. Multiple information points are collected in one sampling cycle. In the figure, f(t) is the original signal to be collected, and f T(t) is the original signal that repeats with a period of T, f △ (t) is T s =T+△t is the periodic pulse trigger signal, △t is the delay step, t p is the effective signal time collected, τ △ is the real-time sampling interval, f s1 (t) is the real-time sampling reconstruction result with △t as the sampling interval, f s4 (t) is the result of fast equivalent sampling reconstruction. Figure 4 It can be seen that unlike the traditional equivalent sampling, the improved equivalent sampling collects K points each time, so the number of required cycles is reduced to M (=N / K), which greatly improves the information acquisition rate of the system.

[0060] In the process of equivalent sampling, the technology for generating precise time steps is very critical. Whether the time step is accurate directly affects the degree of consistency between the reconstructed waveform and the real waveform. In the glacier thickness radar, the time step signal is generated by the timing control module and implemented by a programmable delay chip. The minimum time delay step Δt of the chip is 2ps, and the maximum achievable delay is about 2ns. Δt changes with the operating temperature and the operating voltage. For some key circuits of the timing control module, voltage stabilization and temperature control measures are adopted to monitor the working environment of the programmable delay chip. Based on the PID algorithm in control theory, the working environment parameters are dynamically adjusted using voltage stabilization and temperature control devices, thereby improving the stability and accuracy of the time step.

[0061] Step 4: Real-time positioning of the collected radar wave signal: In the process of receiving the echo, synchronize RTK to obtain the corresponding position coordinates. After the RTK information is collected, the longitude and latitude information therein can be used to realize the dynamic display of the measurement trajectory on the electronic map. At the same time, the longitude and latitude information is stored in the database table corresponding to the measurement data. Combined with GIS (Geographic Information System), the location information corresponding to the measurement data can be displayed;

[0062] RTK positioning is a real-time dynamic positioning technology based on carrier phase observation, which can provide real-time three-dimensional positioning results of mobile stations in a specified coordinate system with centimeter-level accuracy. The RTK receiver is connected to the control computer on the ground through the RS232 serial port or wireless communication interface to collect and store data and control and adjust the flight attitude and corresponding positioning data.

[0063] Step 5, ice thickness data processing: The processing includes the following sub-steps:

[0064] Sub-step 5.1, filtering and gain adjustment, improving radar echo quality and measurement accuracy;

[0065] There are various interference sources in the actual detection environment, such as electromagnetic interference generated by surrounding electrical equipment, surface clutter, etc. If these interference signals are not processed, they will be mixed with the target echo signal, seriously affecting the identification and measurement of underground targets. Filtering is to filter out interference signals outside the target signal frequency band according to the frequency characteristics of the signal. A bandpass filter is used to accurately select the frequency band where the target signal is located to make the useful signal clearer. Since the ground penetrating radar signal will attenuate during the propagation under the ice layer, the echo signal strength will gradually weaken as the detection depth increases. Nonlinear gain adjustment is used to compensate for echo signals at different depths, enhance the strength of deep target echo signals, and make them easier to compare and analyze with shallow target echo signals in terms of strength. Through filtering and gain adjustment, it is ensured that the signals within the entire detection depth range can be accurately processed, thereby improving the measurement accuracy.

[0066] Sub-step 5.2, convert the two-way travel time of the radar signal into the ice thickness H to form the ice thickness layer on the survey line. The calculation formula is:

[0067]

[0068] Where: c is the propagation speed of electromagnetic waves in vacuum (0.29979m / ns); t is the two-way travel time of electromagnetic waves in the measurement medium; ε is the relative dielectric constant of the medium; v is the propagation speed of radar pulses in the medium; d is the distance between the radar transmitting antenna and the receiving antenna;

[0069] Sub-step 5.3, the obtained ice thickness data is marked on the survey line, and the data includes ice thickness information and corresponding geographical location coordinates (such as longitude and latitude) to form a point layer;

[0070] Sub-step 5.4, superimposing the ice thickness point layer to the glacier boundary, superimposing the thickness point layer data on the glacier map according to the geographic coordinates, so that the thickness information on the ice layer measurement line is added to the glacier boundary;

[0071] Sub-step 5.5, interpolate to obtain the ice thickness distribution map of all locations within the glacier boundary; in actual measurement, only the ice thickness at the survey line location can be obtained. Using the linear interpolation method, based on the existing ice thickness data on the survey line, the ice thickness of the area between the survey lines is estimated, thereby obtaining the ice thickness distribution data of all locations within the glacier boundary;

[0072] Step 6, glacier reserve estimation: Estimate the glacier ice reserve based on the obtained glacier body and ice thickness distribution map using block volume method, thickness integral method, etc.

[0073] Based on the obtained glacier body and ice thickness distribution map, the block volume method is used to estimate the glacier ice reserves. The glacier area is divided into multiple area units according to thickness. The volume of each unit is approximately the unit area multiplied by the corresponding glacier thickness, and then the volume of all units is summed. The ice reserve estimation specifically includes the following sub-steps:

[0074] Sub-step 6.1, based on the ice thickness distribution set on the glacier radar survey line, form an ice thickness contour distribution map; divide it into 5m or 10m intervals along the elevation direction, and draw a map in GIS or other software;

[0075] Sub-step 6.2, estimate the bottom area occupied by the glacier body at the same thickness in the GIS map;

[0076] Sub-step 6.3, the volume of the glacier block is obtained by multiplying the bottom area by the partition thickness; the volumes of the blocks at other ice thickness heights can be obtained in turn; and the total volume of the measured glacier body can be obtained by cumulative summation.

[0077] Application examples:

[0078] Dongkemadi Glacier is located in Tanggula Mountain. It is a continental glacier. It is a composite valley glacier formed by the confluence of a south-facing main glacier (Big Dongkemadi Glacier) and a southwest-facing branch glacier (Small Dongkemadi Glacier). The glacier area is 16km 2 , 5.4km long, with the terminal elevation at 5275m. Xiaodongkemadi Glacier (33°04'N, 92°04'E) is a typical continental glacier, with the terminal elevation at 5420m, the highest point at 5919m, and the glacier surface mainly distributed between 5550-5790m above sea level.

[0079] The drone glacier thickness radar system described in this embodiment was used to measure the Xiaodongkemadi Glacier at an altitude of more than 5000m. The drone flew along the survey line for measurement (only for verification and explanation), with a total length of about 4km, and obtained a more detailed distribution of glacier thickness. The measured data show that the echo signal of the ice-rock interface of the glacier radar is clear, and the thickness of the glacier can be accurately obtained. The glacier thickness distribution of each survey line measured on the Xiaodongkemadi Glacier shows an overall trend of larger ice thickness in the middle, gradually decreasing on both sides, and smaller ice thickness at low altitudes. As the altitude increases, the ice thickness gradually increases. When the altitude exceeds 5600m, the rate of increase in glacier thickness gradually decreases. The maximum measured glacier thickness is 125.54m (altitude 5609m).

[0080] The three-dimensional distribution and thickness contour map of Xiaodongkemadi Glacier were obtained by difference processing using the minimum curvature method. Based on this map, the glacier reserves within the measurement range were calculated using the block volume method, and the measured glacier reserves were found to be 69.467 million m 3 .

[0081] Finally, it should be noted that the above is only used to illustrate the technical solution of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the preferred arrangement scheme, a person skilled in the art should understand that the technical solution of the present invention (such as the form of the UAV, the form of the ground penetrating radar, the installation method on the UAV, the sequence of steps, etc.) can be modified or replaced by equivalents without departing from the spirit and scope of the technical solution of the present invention.

Claims

1. A method for estimating glacier ice reserves based on a flying ice radar, the system used in the method comprising: An unmanned aerial vehicle and a glacier thickness measuring radar mounted on the unmanned aerial vehicle, and a control computer arranged on the ground; the glacier thickness measuring radar comprises an ultra-wideband antenna unit, a transmitter, a receiver, an RTK positioning device capable of achieving centimeter-level positioning accuracy, and a signal processing device; the ultra-wideband antenna unit is a linear dipole transmitting antenna and a receiving antenna; the radar pulse signal of the transmitter is generated by using avalanche triodes and step recovery diodes to form a differential structure and a Marx cascade structure to form a carrier-free pulse generating circuit, and adopts a multi-tube joint working mode; the receiver is provided with a sensitivity time control circuit and a time-varying gain amplifier circuit to enhance the reflected signal of the target and increase the detection depth of the radar; The method is characterized in that the steps are as follows: Step 1: Obtaining large-scale glacier topography: Either using public satellite remote sensing data or digital models, or using topographic mapping equipment for on-site mapping, obtain the boundary and surface elevation data of the glacier being studied; Step 2, transmitting radar waves for measuring ice thickness within the boundary of the glacier: planning a survey line within the boundary of the glacier, and setting a flight path in the UAV according to the survey line; starting the UAV to cruise according to the pre-set flight path, turning on the transmitter of the glacier thickness radar, using avalanche triodes and SRD devices to generate carrier-free pulses that meet the requirements of radar detection of glaciers, and using a multi-tube cascade working mode to generate ultra-wideband low phase noise pulse transmission signals, and transmitting detection radar waves through a linear dipole antenna; Step 3, receiving and measuring the radar echo reflected in the glacier: the receiver uses a sensitivity time control circuit and a time-varying gain amplifier circuit to increase the amplifier gain of the receiver in a nonlinear form as the pulse propagation time increases, so as to enhance the reflected signal of the target and increase the detection depth of the radar; the processing of the received signal adopts an equivalent sampling method to reconstruct the radar echo signal; Step 4: Real-time positioning of the collected radar wave signal: In the process of receiving the echo, synchronize RTK to obtain the corresponding position coordinates. After the RTK information is collected, the longitude and latitude information can be used to track and dynamically display the measurement trajectory on the electronic map. At the same time, the longitude and latitude information is stored in the database table corresponding to the measurement data. Combined with GIS, the location information corresponding to the measurement data can be displayed; Step 5, ice thickness data processing: The processing includes the following sub-steps: Sub-step 5.1, filtering and gain adjustment, improving radar echo quality and measurement accuracy; Sub-step 5.2, convert the two-way travel time of the radar signal into the ice thickness H to form the ice thickness layer on the survey line. The calculation formula is: Where: c is the propagation speed of electromagnetic waves in vacuum; t is the two-way travel time of electromagnetic waves in the measurement medium; ε is the relative dielectric constant of the medium; v is the propagation speed of radar pulses in the medium; d is the distance between the radar transmitting antenna and the receiving antenna; Sub-step 5.3, the obtained ice thickness data is marked on the survey line, and the data includes ice thickness information and corresponding geographical location coordinates to form a point layer; Sub-step 5.4, superimposing the ice thickness point layer to the glacier boundary, superimposing the thickness point layer data on the glacier map according to the geographic coordinates, so that the thickness information on the ice layer measurement line is added to the glacier boundary; Sub-step 5.5, interpolate to obtain the ice thickness distribution map of all locations within the glacier boundary; in actual measurement, only the ice thickness at the measuring line position can be obtained; use the linear interpolation method to estimate the ice thickness of the area between the measuring lines based on the existing ice thickness data on the measuring lines, so as to obtain the ice thickness distribution data of all locations within the glacier boundary; Step 6, glacier reserve estimation: The reserve estimation includes the following sub-steps: Sub-step 6.1, based on the ice thickness distribution data on the glacier radar survey line, form an ice thickness contour distribution map; divide it into 5m or 10m intervals along the elevation direction, and draw a map in GIS; Sub-step 6.2, estimate the bottom area occupied by the glacier body at the same thickness in the GIS map; Sub-step 6.3, the volume of the glacier block is obtained by multiplying the bottom area by the partition thickness; the volumes of the blocks at other ice thickness heights can be obtained in turn; and the total volume of the measured glacier body can be obtained by cumulative summation.

Images