Nothing Special   »   [go: up one dir, main page]
Next Article in Journal
Landslide Mapping and Monitoring Using Persistent Scatterer Interferometry (PSI) Technique in the French Alps
Next Article in Special Issue
SCoBi Multilayer: A Signals of Opportunity Reflectometry Model for Multilayer Dielectric Reflections
Previous Article in Journal
Detecting Change in Forest Structure with Simulated GEDI Lidar Waveforms: A Case Study of the Hemlock Woolly Adelgid (HWA; Adelges tsugae) Infestation
Previous Article in Special Issue
A Modified Model for Electromagnetic Scattering of Sea Surface Covered with Crest Foam and Static Foam
You seem to have javascript disabled. Please note that many of the page functionalities won't work as expected without javascript enabled.
 
 
Search for Articles:
Title / Keyword
Author / Affiliation / Email
Journal
Article Type
 
 
Advanced Search
 
Section
Special Issue
Volume
Issue
Number
Page
 
Journals
Remote Sensing
Volume 12
Issue 8
10.3390/rs12081300
remotesensing-logo
Submit to this Journal Review for this Journal Propose a Special Issue
► Article Menu

Article Menu

share Share announcement Help format_quote Cite question_answer Discuss in SciProfiles thumb_up
...
Endorse textsms
...
Comment
first_page
settings
Order Article Reprints
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Mapping and Assessment of Tree Roots Using Ground Penetrating Radar with Low-Cost GPS

by
Lilong Zou
1,*,
Yan Wang
2,
Iraklis Giannakis
1,
Fabio Tosti
1,
Amir M. Alani
1 and
Motoyuki Sato
3
1
School of Computing and Engineering, University of West London (UWL), London W5 5RF, UK
2
Department of Earth and Space Sciences, Southern University of Science and Technology, Shenzhen 518055, China
3
Center for Northeast Asian Studies, Tohoku University, Sendai 9808576, Japan
*
Author to whom correspondence should be addressed.
Remote Sens. 2020, 12(8), 1300; https://doi.org/10.3390/rs12081300
Submission received: 17 March 2020 / Revised: 14 April 2020 / Accepted: 18 April 2020 / Published: 20 April 2020
(This article belongs to the Special Issue Electromagnetic Modeling in Microwave Remote Sensing)
Browse Figures
Graphical abstract
"> Figure 1
<p>(<b>a</b>) Surveying scenario of the Metasequoia trees using the proposed 3D GPR system. (<b>b</b>) The low-cost GlobalSat GPS receiver used in this investigation.</p> "> Figure 2
<p>Time-domain reflectometry (TDR) measurements of the soil moisture for the estimation of the subsurface velocity in the site.</p> "> Figure 3
<p>Raw radargram of the tree root investigation using a 500 MHz RAMAC shielded antenna system.</p> "> Figure 4
<p>(<b>a</b>) The location of the GEONET stations in Japan; (<b>b</b>) location of the investigation site and nearest GEONET Station.</p> "> Figure 5
<p>(<b>a</b>) GPS receiver geometry with a reference station; (<b>b</b>) flowchart of the GPS receiver bias estimation and removal by use of the Kalman filter.</p> "> Figure 6
<p>Moving trajectories of antennas recorded by the GlobalSat GPS receiver. (<b>a</b>) The record by the GlobalSat GPS receiver; (<b>b</b>) post-processing results.</p> "> Figure 7
<p>Migrated vertical profiles along the survey direction. (<b>a</b>) Migrated profile at 0.8 m cross-survey direction; (<b>b</b>) migrated profile at 1.6 m cross-survey direction; (<b>c</b>) migrated profile at 2.4 m cross-survey direction; (<b>d</b>) migrated profile at 3.2 m cross-survey direction.</p> "> Figure 8
<p>Migrated horizontal slices at different depths. Red lines indicate the detected roots in the migrated data set. (<b>a</b>) Migrated slice at 16-cm depth; (<b>b</b>) migrated slice at 19-cm depth; (<b>c</b>) migrated slice at 27.5-cm depth.</p> "> Figure 9
<p>Migrated horizontal slices at different depths. Red lines indicate the detected roots in the migrated data set. (<b>a</b>) Migrated slice at 42-cm depth; (<b>b</b>) migrated slice at 59-cm depth; (<b>c</b>) migrated slice at 79.5-cm depth.</p> "> Figure 10
<p>The coarse root detection in a 3D migrated cubic: (<b>a</b>) B-scan profile which contains a tree root. (<b>b</b>) Root extension tracking with the −3 dB positions of a local peak (survey direction). (<b>c</b>) Root extension tracking with the −3 dB positions of a local peak (cross-survey direction).</p> "> Figure 11
<p>The reconstructed 3D root system: (<b>a</b>) View from the starting point; (<b>b</b>) view from the ending point.</p> ">
Versions Notes

Abstract

:
In this paper, we have presented a methodology combining ground penetrating radar (GPR) and a low-cost GPS receiver for three-dimensional detection of tree roots. This research aims to provide an effective and affordable testing tool to assess the root system of a number of trees. For this purpose, a low-cost GPS receiver was used, which recorded the approximate position of each GPR track, collected with a 500 MHz RAMAC shielded antenna. A dedicated post-processing methodology based on the precise position of the satellite data, satellite clock offsets data, and a local reference Global Navigation Satellite System (GNSS) Earth Observation Network System (GEONET) Station close to the survey site was developed. Firstly, the positioning information of local GEONET stations was used to filter out the errors caused by satellite position error, satellite clock offset, and ionosphere. In addition, the advanced Kalman filter was designed to minimise receiver offset and the multipath error, in order to obtain a high precision position of each GPR track. Kirchhoff migration considering near-field effect was used to identify the three-dimensional distribution of the root. In a later stage, a novel processing scheme was used to detect and clearly map the coarse roots of the investigated tree. A successful case study is proposed, which supports the following premise: the current scheme is an affordable and accurate mapping method of the root system architecture.

Graphical Abstract">

Graphical Abstract

1. Introduction

Environmental issues such as the conservation of natural heritage and ancient trees have become priority objectives of urgent protection [1,2]. Very limited knowledge is present on how the elements causing the rapid death of entire forests interact with each other. The quality and distribution of roots can be used as reliable diagnostic parameters for early tree decline [3,4,5]. Therefore, if the depth of the root and its extension to the surrounding area can be estimated, more effective preconditions can be provided to guide a wide range of research decisions [6].
Destructive methods can provide accurate mapping of the root systems. Nonetheless, intrusive approaches are often impossible to implement in the field and time-consuming. In addition, digging and trenching can affect the surrounding of trees and may cause irreversible damage [7]. Destructive methods are limited in space and can only provide local information on the tree conditions. Instead, non-destructive testing (NDT) methods are becoming more and more common in forestry and arboriculture applications due to the high productivity and reliability of the information provided. Based on the aforementioned factors, GPR technology has been proven to be one of the most effective NDT tools. GPR is characterized by a high versatility, a fast data collection, and the provision of reliable results, at relatively limited costs. GPR has been widely used in various disciplines, such as pavement analysis [8,9], archaeological investigations [10], mine detection [11] and civil and environmental engineering applications [12]. The application of GPR in forestry sciences is usually related to tree trunk assessment, root mapping, and the evaluation of the soil-tree interaction. The first application of GPR in tree root research can be traced back to 1999 [13], referring to the mapping of tree root systems. Since then, research has focused on the evaluation of the root diameter [14,15,16] in urban areas, the functions of the root system architecture such as biomass [17,18,19,20], roughness [21,22,23] and root mode [24,25,26].
In recent years, with the further development of GPR, especially the combination with advanced positioning systems, it was possible to detect the 3D structure in the subsurface with a 3D full-resolution, such as in the case of the tree root system presented in [27]. In [27], a 3D GPR system coupled with a high precision position scheme was used for mapping tree roots. This high-precision framework is based on three indoor GPS transmitters and uses a laser receiver fixed on the GPR to track the location of the each GPR trace. At the cost of expensive equipment and time-consuming field measurement, it can achieve the accuracy of hundreds of microns. In [28], the rotary laser positioning technology has been applied in 3D GPR measurements. The positioning system enables a centimeter accuracy to be achieved by using small detectors attached to the GPR antennas. In [29], another high-precision positioning system, self-tracking total station, was used for 3D GPR field measurements purposes. The self-tracking total station needs a cable connection, so it is not a flexible in-field application. Although the wireless communication between GPR and self-tracking total station (TTS) system was developed to avoid cable connection, the crosstalk effect and the time synchronization required between the GPR and the positioning system have a significant impact on the data quality of GPR. Other instruments such as charge-coupled device (CCD) cameras, real-time kinematic GPS (RTK GPS) and gyroscope can also be used for actual GPR measurements. However, due to the high costs and the large demand for calculation and operation, the applicability of these methods in real-life applications is limited.
Within this context, the aim of this study was to develop a methodology to process the data acquired by a cost-effective 3D GPR system. Cost-effectiveness is intended in a way that the system has a low system complexity and a higher efficiency. The methodology is coupled with a set of optimisation algorithms in order to achieve fast and affordable detection tools fine-tuned for the detection of roots.
The rest of the paper is divided into five main sections. Section 2 describes the aims and objectives of this paper. Section 3 introduces the survey site and the data acquisition in the field. This section also discusses the low-cost GlobalSat GPS Receiver and the B-scan outputs collected in the measurements. Section 4 presents the data processing methodology that includes the pre-processing and the post-processing of the tracking position, and the Kirchhoff migration. Section 5 discusses a novel approach to isolate and highlight the coarse roots of the tree, providing a map about their 3D distribution. Section 6 summarises the main findings and conclusions of this research.

2. Aims and Objectives

The main aim of the current research was to provide an effective, economical, and full-resolution mapping for the assessment of the roots system.
To that extent, the following objectives were identified: First, to develop an algorithm that could provide high-precision position information for each GPR trace by a low-cost GPS positioning system. Second, to develop a novel data-processing approach for achieving 3D full-resolution visualisation of shallow and deep region tree roots.

3. Survey Site and Data Acquisition

The tree roots field measurements were carried out over a large open lawn in front of the Tohoku University Centennial Hall–Kawauchi Hagi Hall. The Hall is located in downtown Sendai (Japan) among other international venues and cultural facilities. The site was once the place where the U.S. military was stationed after World War II. After that, the Kawauchi Hagi Hall and the surrounding facilities were built in this area. The survey site was located in a rectangular area measuring 10 m × 4 m in the large open lawn, about 6-m beside the surrounding trees. The latter are typical Metasequoia glyptostroboides, with an average diameter at breast height of around 50 cm and a tree height that can approach 20 m. Figure 1a shows the scenario of the investigated site. Metasequoia glyptostroboides belong to a shallow-rooted species. The main roots of trees are underdeveloped, whereas the lateral roots or adventitious roots grown by radiation are much longer. Most of the roots are distributed on the surface of the soil. Lateral roots are mainly distributed between 0.2–1 m, with some tree roots being exposed to the ground. The white dashed line with an arrow in Figure 1a indicates the investigation direction crossing the root distribution that was used for the data collection stage.
For the purpose of this investigation, a 500 MHz RAMAC shielded bowtie antenna was employed. Frequency characteristics of the antenna is summarised in Table 1. Using a low-cost GlobalSat GPS receiver (Figure 1b) in combination with a GPR system, the positional coordinates of the moving antenna was able to be recorded in the collected GPR track. In order to obtain the exact location of the relevant GPR trace, the local GEONET station and the GNSS information were needed for post-processing of the position data collected by the low-cost GPR on the site.
In this study, to make the maximum amount of clear three-dimensional views of the tree roots and understand how the root distributed over the subsurface, high-density GPR data acquisitions are usually required. According to [27,30], the sampling rate Δ x must be equal to or less than the Nyquist sampling rate Δ x s along the direction of investigation and cross-investigation for the full-resolution 3D imaging:
Δ x ≤ Δ x s = λ s i n θ
where θ is the main lobe angle of the antenna that was used in the investigation and λ is the wavelength of the subsurface material. Since the main lobe angle of the GPR antenna used in the near field is usually greater than 60°, the maximum space interval should be equal to or less than one-quarter wavelength of the antenna bandwidth on the subsurface [27].
In order to evaluate the dielectric constant of the soil, the time-domain reflectometry (TDR) technique was used at several locations in the measurement area, as shown in Figure 2. The average value of the volumetric water content in the subsurface material was 42%, so that the subsurface could be classified as a wet soil. The relative permittivity measured near the soil surface was 27 which was calculated using the equation described in [31,32,33]. The corresponding subsurface velocity was 0.06 m/ns.
According to the parameters of the 500 MHz RAMAC shielded antenna in Table 1, the underground velocity v is equal to 0.06 m/ns, and the spatial sampling interval is equal to or less than 5.5 cm. Therefore, the interval along the investigation direction should set less than 5.5 cm to meet the conditions of full-resolution 3D imaging. Under these circumstances, the system could provide a 3D subsurface image with 5.5 cm of horizontal resolution and 6.7 cm of vertical resolution in this survey scenario. Moreover, it is worth mentioning that an excess of 900 m of surveys with the GlobalSat GPS Receiver were performed in approximately half an hour. In particular, more than 90 surveys have been carried out across the 4-m cross survey direction. The raw radargram of the whole tree root survey is shown in Figure 3. A superposition of multiple responses could be seen up to 30 ns.

4. Data Processing and Methodology

The processing presented in this study includes four sequential stages. In the pre-processing stage, time correction and a signal noise filter were applied to increase the clutter rate of the entire data. Then, by using the local GEONET station, Ultra-Rapid ephemeris, precise clock data, and an optimisation algorithm were able to obtain accurate location information, which was used. Next, the near-field Kirchhoff migration was used to map the distribution of the coarse roots. In the last stage, a novel processing scheme was used to identify and clearly map the coarse roots of the investigated tree. All the algorithms in each stage were implemented in MATLAB.

4.1. Radargram Signal Processing

Prior to any interpretation attempt, the raw data were subject to a pre-processing step in order to reduce undesired reflections and increase the overall signal-to-clutter ratio. In particular, the pre-processing step consisted of:
(a) Time-zero adjustment
This step is applied to adjust the initial position of the surface reflection in GPR signals (the time when the radar signal spreads out from the antenna and enters the ground is known as “zero time”). In order to move all traces to the zero-time position before other processing methods, time-zero adjustment must be applied.
(b) Background removal
Cross-coupling between the transmitter and the receiver as well as a ringing noise and multiple reflections can mask the less dominant reflections from the root. In order to mitigate such unwanted signals, background removal methods should be applied. Subtracting the background signal proved very effective for the investigated case study.
(c) Band-pass filtering
In order to improve the quality of the data, a band-pass filter is applied to eliminate different types of noise in this stage. In addition, filters are also used to extract useful information and hidden patterns from the recorded data. In order to remove the direct current offset and suppress the high-frequency noise, a band-pass filter was employed with the cut off frequencies F L o w and F H i g h indicated in Table 1.

4.2. Post-Processing of the Tracking Position

After pre-processing the GPR raw data, the post-processing will be performed on the GPR trace position data collected by the low-cost GPS receiver. So that the high precision position of each GPR trace can be applied for the 3D full resolution imaging for further investigation. GPS provides continuous global position, time and navigation services. The satellites send navigation and observation data to the receivers [34,35,36]. As the GPS receiver is affected by a number of different parameters, estimating the receiver position is not a trivial task. These parameters include GPS receiver clock bias, troposphere, ionosphere, multipath propagation and so on [37,38].
Each of these factors can be modelled separately. The basic pseudo-range equation is given as [39]:
P m = ρ g + c ( Δ t r − Δ t s ) + I ε + T ε + ε m u l ,
where c is the speed of light, T ε is the error produced by the troposphere, I ε is the error produced by the ionosphere and ε m u l is the multipath error. Δ t s is the offset generated by the satellite bias and Δ t r is the receiver clock bias. ρ g is the geometric distance from the GPS satellite to the receiver which can be given by:
ρ g = ( x r − x s ) 2 + ( y r − y s ) 2 + ( z r − z s ) 2 ,
where ( x r , y r , z r ) and ( x s , y s , z s ) are the positional coordinates of the receiver and a given satellite, respectively.
In Equation (2), the error term of correction is divided into two groups. The troposphere error, the ionosphere error and the satellite clock bias are set into one group, represented as R s . The receiver clock bias and the multipath error are set into another group, represented as R r . Therefore, Equation (2) can be formulated as:
P m = ( x r − x s ) 2 + ( y r − y s ) 2 + ( z r − z s ) 2 + R s + R r ,
In Equation (4), the receiver coordinates can be accurately calculated after correctly estimating the items R s and R r . By using signals from at least four satellites at the same time, the GlobalSat GPS receiver can provide real-time single point positioning [40]. The weakness of GPS comes from the satellite position, the satellite clock offset and the ionosphere influence, amongst others [41,42,43]. Nevertheless, the accuracy can be improved by post-processing the received signal.
Three to nine hours after the observation, the observed Ultra-Rapid ephemeris and clock data can be obtained online. This information is released four times a day at 03:00, 09:00, 15:00 and 21:00 UTC. The predicted Ultra-Rapid ephemeris and clock data are also released at the same time, which can be obtained in advance. Instead, information from Rapid 1 takes longer (from 17 to 41 h) to get online. It is released every day at 17:00 UTC [44,45,46,47,48].
International GNSS Service (IGS) is a cooperation with many organisations and countries. This collects data from more than 300 continuously operating reference stations around the world. Based on these data, it develops precise satellite ephemeris and clock solutions. The processing phase involves up to eight IGS analysis centers and the results are freely distributed by IGS. To this effect, this service allows researchers to post-process the observations based on the information provided by the IGS. Compared with a broadcast ephemeris and clock correction, the aforementioned approach is more accurate since these data reflect the precise position and clock offset of the satellite during the actual measurement. These data are classified into several categories and discussed in the literature [49,50,51].
In general, differential processing techniques depend on at least two receivers standing at control stations (known as a “bases”) with a known location. In Japan, there exist around 1200 GEONET stations that record data every 30 s, as shown in Figure 4a. Therefore, it is possible to know the magnitude of the basic position error of the receiver. If the difference between the deviation of the reference station and the deviation of the flow station can be found, the position error at the reference line can be calculated. With the differential process, a correction is generated. The procedure can reduce the position error of unknown points. In general, this method can provide the location with sub-meter resolution from single-frequency pseudo-range observations. However, this accuracy level still cannot meet the requirements of current applications.
It is also worth noting that this application does not require a global solution. Only the 2D relative position of each GPR acquisition point needs to be estimated accurately. For this purpose, we selected a reference GEONET station (as shown in Figure 4b) 14 km away from the station. Assuming that the troposphere error T ε , the ionosphere error I ε and the satellite clock bias of the reference station were the same as those of the experiment site, the R s terms in Equation (4) can be directly estimated. Moreover, the receiver bias ρ g and the multipath error ε m u l can be estimated and suppressed by designing an advanced extended Kalman filter.
In Equation (4), R s is mainly caused by the substantial errors of the positions of the satellites, the clock offset, the ionosphere correction etc. Using the record data from the reference station, a correction matrix was designed with the accuracy of the satellite’s positions, clock offset and ionosphere information, which can be downloaded afterwards from the International GNSS Service (IGS) database in order to minimise R r e s i d u a l . The measurement equation could be written as:
P r = H 0 [ x r e f y r e f z r e f ] + R r e s i d u a l ,
where P r is the post-processed position of the reference station, x r e f ,   y r e f ,     z r e f are the recorded coordinates of the data at the reference station, H 0 is the design matrix that contains the satellite and the reference station receiver clock biases plus the ionospheric and tropospheric effects and another minor correction term. Lastly, R r e s i d u a l is the residual error.
Based on the aforementioned, the least square term of the residual error R s in the experiment period can be written as:
R L S = R r e s i d u a l T R r e s i d u a l .
Consequently, the optimally-designed matrix H 0 can be calculated as:
min R L S → 0 ‖ H 0 T H 0 x − H 0 T y ‖ 2 = 0
and Equation (4) can be written as:
P m = H [ x r y r z r ] + R s ≅ H 0 [ x r y r z r ] + R s ,
where P m is the post-processed position of each acquisition point, x is the recorded data at each GPR acquisition point, H is the design matrix, which can be replaced by H 0 (as described in Figure 5a) and R s is the receiver clock bias and other minor correction terms such as the multipath error and the residual error.
The R s term is estimated and removed based on the conventional Kalman forward-backward filter. The estimated offset of the receiver clock located at the tracking station does not represent the actual offset due to the fact that the obtained observation data were pre-processed [52,53,54]. An apriority bias of the GPS receiver was used in the beginning to estimate the term R s . Finally, all measurements were then corrected by the terms H and R s .
In order to update the filter over time, the predictive state vector of the system model should be used. Since the GPS receiver bias is assumed to be a linear function of time, the drift and all the other parameters can be considered as constants. It is worth noting that the drift is not strictly the same, as it will change slowly with time. Figure 5b depicts the complete flowchart of the proposed methodology consisting of six distinct and sequential steps:
(a) Forward Filter Initialisation
The coarse value of the GlobalSat GPS receiver is used as a prior value for the bias and the drift, with all the other elements set to zero. In addition, the noise of the filter needs to be set in this step.
(b) GPS Bias and Multi-Pass Error Constraint
In this step, the GPS bias and multi-pass error will propagate towards the current filter state.
(c) Kalman Filter Measurement Update
After the filter run is finished, a state vector computes the mean of the forward and backward results of the filter state.
(d) Results and Covariance Matrix
In this step, the results are stored to be used in the subsequent filter at this stage. In the meantime, the covariance matrix is needed to be evaluated. Then, the smoothing factor should be calculated.
(e) Judgement
The procedure is repeated until the processing reaches the end time, or until it meets the constraint condition. Here, Allen variance σ a is used to obtain a rough estimate of the expected receiver bias and residual error. The Allan variance σ a is a measure of frequency stability in clocks, oscillators and amplifiers [37].
(f) Smooth operator
After all the procedures are finished, the Kalman filter parameters are stored and used to design the smooth operator.
Figure 6a shows the moving trajectories of the antennas recorded by the GlobalSat GPS receiver. During the investigation, the effects of the ionosphere had little to negligible effects to the overall precision of the algorithm. In fact, the drift of the GPS sensor itself was proven to have the biggest impact on the overall accuracy. The results after applying the Kalman filter and after converting the data to a local Cartesian coordinate system are shown in Figure 6b. The shape of the trajectory was consistent with the actual GPR system movement.

4.3. Kirchhoff Migration

In order to detect the three-dimensional distribution of the tree roots, a high-resolution and accurate imaging algorithm was implemented in this section. The underground 3D geometry can be reconstructed from the results of dense measurements. Three-dimensional migration is required to focus reflection and diffraction and collapse the recorded hyperbolas to their origin [55,56].
In this research, we applied the Kirchhoff migration, which is based on the generalised wave-equation [57,58]. Mathematically, the Kirchhoff migration method can be formulated in the time domain as follows:
U ( r , t ) = 1 8 π 2 ∫ S 0 ∂ R ∂ n [ 1 v R ∂ ∂ t U ( r 0 , t − τ ) + 1 R 2 U ( r 0 , t − τ ) ] d S 0 .
where S 0 represents a source-free surface; n is the unit vector normal to S 0 ; r is the positional vector of the estimated wave field; r 0 is the positional vector of the integration point; v is the propagation speed in media, U represents the wave field within the surface S 0 and
R = | r − r 0 | ,
τ = R v .
Normally, the far-field approximation ω / R v ≫ 1 / R 2 is applied in (9). Therefore, the second term of Equation (9) can be omitted. Since we were oriented in the near-field target, such an assumption was no longer suitable for the purpose of this research. To mitigate that, the entirety of Equation (9) was used in the current study.
As described in Section 3, TDR measurements were conducted at different points within the surveying area. The average of the results turned out to provide a wet soil with a volumetric water content of 42% found using the Topp equation. The measured relative permittivity in the near-surface of the soil was 27 with a corresponding subsurface velocity of 0.06 m/ns. Since most of the roots were located in the shallow region of the subsurface, a constant velocity of the medium was assumed in this paper.

5. Results and Discussion

The previous section presented the processing methodology of our proposed strategy for tree root mapping and assessment. In this section, the results concerning the 3D distribution of the roots system from the field measurements will be discussed.

5.1. 3D Migration Result

The 3D migration results were obtained by applying the Kirchhoff migration discussed above. One-twentieth of the resolution was applied as the grid size during the migration progressing. Figure 7 shows the migrated vertical profiles along the survey direction at different cross-survey directions (0.8 m, 1.6 m, 2.4 m and 3.2 m). Most of the reflections are located in the upper region from 0.2 m to 0.5 m. In Figure 7c,d, there is a deep reflection located around 0.7 m depth between the section 9–10 m along the survey direction. Figure 8 and Figure 9 show the migrated horizontal slices at different depths (0.16 m, 0.19 m, 0.275 m, 0.42 m, 0.59 m and 0.795 m). By displaying animation frames of slices with different migration depths, the root reflections at each horizontal slice are clearly visible.
Based on the migration results, it was also obvious that most of the reflections from the roots were located in the upper region, as shown in Figure 8. In the figures, some areas with a strong amplitude showed the bright blocks feature. These bright blocks appear to be randomly distributed, and they can be mainly attributed to the heterogeneity of the soils, rocks and similar features in the subsurface. Hence, it is difficult to discriminate coarse roots by using GPR from these output types.
There are many reasons that limits the applicability and the accuracy of GPR for detecting tree roots. For example, if a cluster of fine roots surrounds a deeper coarse root, it is difficult to detect the coarse root. Moreover, if the distance between two coarse roots is less than the horizontal and vertical resolution of the observation system at the same time, then two coarse roots will be considered as a unique root. In addition, roots that are parallel to the polarisation of the antenna will result in weak reflections that will be masked by unwanted clutter and noise. Based on the aforementioned factors, the root can be detected according to the reflection intensity, direction and surrounding environment from the migrated slices. The distribution of tree roots must be mapped on a three-dimensional environment, as the reflection amplitude in the same depth slice varies along the root. Therefore, a coherent 3D approach is necessary in order to effectively locate and track the distribution of the roots.

5.2. 3D Roots Detection Result

The root distribution can be tracked from the image index. The coarse tree roots normally can be displayed continuously in the 3D migrated data. Therefore, these image indexes can be associated with the root. In this paper, high-energy regions inside a 3D cube above a certain threshold were used to extract root locations. In Figure 10a, a B-scan profile which contains a tree root is shown. The profile is perpendicular to the root extension direction. In this migrated cube, a 1.5 m × 1.5 m × 0.3 m (survey direction × cross-survey direction × depth) block was selected. The peak position of the cube was found and a vertical curve (depth) and two horizontal curves (survey and cross-survey directions) of the peak were used for next step analysis. Additionally, the −3 dB width of the peak is defined as the root in each direction, as it is shown in Figure 10b,c.
Overall, the accuracy of the root detection using pixel intensities in the image was slightly worse than extracting the feature directly from the waveform. The –3-dB threshold did not remove all the weak reflections, such as those from the heterogeneous soils, rocks, and fine roots. To avoid false detection results, future processing should also include discarding those targets that do not have continuous distribution.
The reconstructed 3D root system is illustrated from two different viewpoints in Figure 11. Six coarse roots across the survey direction can be clearly detected. The false detection results were discarded. Root 1, 2, 3 and 4 were located at the 30-cm region below the subsurface and their distributions were confirmed with a thin measuring flower pole, which can be inserted into the soil and record the depth of these targets. Root 5 and 6 were located at 0.5 m and 0.7 m, respectively. In the field data set, it was impossible to detect all the coarse roots. The characteristics of tree roots radiating along the trunk provided an assumption for identifying the extension of the roots. On the contrary, fine roots always gather on heterogeneous soil and stone. This is because it is not conducive to the growth of trachyte in these areas. Since there was no marked difference between the fine roots and the surrounding environment, it was difficult to find these targets.

6. Conclusions

In this paper, a low-cost GPS receiver combined with a 500 MHz RAMAC shielded antenna was used to obtain a full-resolution 3D image of the root system of several larch trees. By using the low-cost GPS receiver and advanced processing algorithms, an effective and affordable way for mapping tree root systems was proposed. This equipment combination was used to carry out a full-resolution survey over a 10 m × 4 m area in approximately half an hour, which proved the time efficiency of the approach proposed. In addition, a dedicated post-processing methodology was implemented to enhance the accuracy of each GPR trace recorded by the low-cost GPS receiver. After that, a Kirchhoff migration approach considering the near field effects was applied to accurately map the investigated root system. Lastly, the 3D distribution of coarse roots was mapped by considering the local magnitude in a 3D migrated cube. As a result, six large coarse roots were successfully identified. Although it is not possible to have full detection of the entire root system (i.e., coarse and fine roots), this research demonstrates the potential of combining low-cost GPS and GPR antenna systems to achieve a full-resolution 3D GPR survey in a time-effective and affordable manner.

Author Contributions

L.Z. analysed the data, proposed the algorithm and wrote the paper. Y.W. contributed to data analysis and provided useful suggestions. I.G., F.T. and A.M.A. contributed in structuring the focus of the paper and the presentation of the results, as well as in editing the language. M.S. contributed much to the experimental design, data collection, and language correction. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Stokes, A.; Fitter, A.H.; Courts, M.P. Responses of young trees to wind and shading: Effects on root architecture. J. Exp. Bot. 1995, 46, 1139–1146. [Google Scholar] [CrossRef]
  2. Coutts, M.P. Root architecture and tree stability. Plant Soil 1983, 71, 171–188. [Google Scholar] [CrossRef]
  3. Habermehl, A. A new non-destructive method for determining internal wood condition and decay in living trees. Part 1. Principles, method, and apparatus. Arboric. J. 1982, 6, 1–8. [Google Scholar] [CrossRef]
  4. Habermehl, A. A new non-destructive method for determining internal wood condition and decay in living trees. II: Results and further developments. Arboric. J. 1982, 6, 121–130. [Google Scholar] [CrossRef]
  5. Daniels, D.J. Ground Penetrating Radar, 2nd ed.; Institution of Engineering and Technology: London, UK, 2004. [Google Scholar]
  6. Harry, M.J. Ground Penetrating Radar: Theory and Application, 1st ed.; Elsevier Science: Amsterdam, The Netherlands, 2009. [Google Scholar]
  7. Vore, S.L.D. Ground-penetrating radar: An introduction for archaeologists. Geoarchaeology 1997, 54, 527–528. [Google Scholar] [CrossRef]
  8. Kaur, P.; Dana, K.J.; Romero, F.A.; Gucunski, N. Automated GPR rebar analysis for robotic bridge deck evaluation. IEEE Control Syst. Lett. 2015, 46, 2265–2276. [Google Scholar] [CrossRef]
  9. Zou, L.; Yi, L.; Sato, M. On the Use of Lateral Wave for the Interlayer Debonding Detecting in an Asphalt Airport Pavement Using a Multistatic GPR System. IEEE Trans. Geosci. Remote Sens. 2020, 58, 1–10, (early access). [Google Scholar] [CrossRef] [Green Version]
  10. Zhou, H.; Sato, M. Archaeological investigation in Sendai Castle using ground-penetrating radar. Archaeol. Prospect. 2001, 8, 1–11. [Google Scholar] [CrossRef]
  11. Sato, M.; Hamada, Y.; Feng, X.; Kong, F.N.; Zeng, Z.; Fang, G. GPR using an array antenna for landmine detection. Near Surf. Geophys. 2004, 2, 7–13. [Google Scholar] [CrossRef]
  12. Slob, E.; Sato, M.; Olhoeft, G. Surface and borehole ground-penetrating-radar developments. Geophysics 2010, 75, 75A103–75A120. [Google Scholar] [CrossRef] [Green Version]
  13. Hruska, J.; Čermák, J.; Šustek, S. Mapping tree root systems with ground-penetrating radar. Tree Physiol. 1999, 19, 125–130. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  14. Butnor, J.R.; Doolittle, J.A.; Kress, L.; Cohen, S.; Johnsen, K.H. Use of ground-penetrating radar to study tree roots in the southeastern United States. Tree Physiol. 2001, 21, 1269–1278. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  15. Barton, C.V.; Montagu, K.D. Detection of tree roots and determination of root diameters by ground penetrating radar under optimal conditions. Tree Physiol. 2004, 24, 1323–1331. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  16. Butnor, J.R.; Doolittle, J.A.; Johnsen, K.H.; Samuelson, L.; Stokes, T.; Kress, L. Utility of ground-penetrating radar as a root biomass survey tool in forest systems. J. Soil Sci. 2003, 67, 1607–1615. [Google Scholar] [CrossRef] [Green Version]
  17. Hirano, Y.; Dannoura, M.; Aono, K.; Igarashi, T.; Ishii, M.; Yamase, K.; Makita, N.; Kanazawa, Y. Limiting factors in the detection of tree roots using ground-penetrating radar. Plant Soil. 2009, 319, 15. [Google Scholar] [CrossRef]
  18. Zenone, T.; Morelli, G.; Teobaldelli, M.; Fischanger, F.; Matteucci, M.; Sordini, M.; Armani, A.; Ferrè, C.; Chiti, T.; Seufert, G. Preliminary use of ground-penetrating radar and electrical resistivity tomography to study tree roots in pine forests and poplar plantations. Plant Physiol. 2008, 35, 1047–1058. [Google Scholar] [CrossRef]
  19. Liu, H.; Huang, X.Y.; Han, F.; Cui, J.; Spencer, B.F.; Xie, X.Y. Hybrid Polarimetric GPR Calibration and Elongated Object Orientation Estimation. IEEE J. Sel. Top. Appl. Earth. Obs. Remote Sens. 2019, 12, 2080–2087. [Google Scholar] [CrossRef]
  20. Cui, X.; Guo, L.; Chen, J.; Chen, X.; Zhu, X. Estimating tree-root biomass in different depths using ground-penetrating radar: Evidence from a controlled experiment. IEEE Trans. Geosci. Remote Sens. 2012, 51, 3410–3423. [Google Scholar] [CrossRef]
  21. Raz-Yaseef, N.; Koteen, L.; Baldocchi, D.D. Coarse root distribution of a semi-arid oak savanna estimated with ground penetrating radar. J. Geophys. Res. Earth Surf. 2013, 118, 135–147. [Google Scholar] [CrossRef] [Green Version]
  22. Li, W.; Cui, X.; Guo, L.; Chen, J.; Chen, X.; Cao, X. Tree Root Automatic Recognition in Ground Penetrating Radar Profiles Based on Randomized Hough Transform. Remote Sens. 2016, 8, 430. [Google Scholar] [CrossRef] [Green Version]
  23. Leucci, G. The use of three geophysical methods for 3D images of total root volume of soil in urban environments. Explor. Geophys. 2010, 41, 268–278. [Google Scholar] [CrossRef]
  24. Tanikawa, T.; Hirano, Y.; Dannoura, M.; Yamase, K.; Aono, K.; Ishii, M.; Igarashi, T.; Ikeno, H.; Kanazawa, Y. Root orientation can affect detection accuracy of ground-penetrating radar. Plant Soil. 2013, 373, 317–327. [Google Scholar] [CrossRef]
  25. Alani, A.M.; Lantini, L. Recent advances in tree root mapping and assessment using non-destructive testing methods: A focus on ground penetrating radar. Surv. Geophys. 2020, 41, 605–646. [Google Scholar] [CrossRef]
  26. Tosti, F.; Bianchini Ciampoli, L.; Brancadoro, M.G.; Alani, A. GPR applications in mapping the subsurface root system of street trees with road safety-critical implications. Adv. Trans. Stud. 2018, 44, 107–118. [Google Scholar]
  27. Zhu, S.; Huang, C.; Su, Y.; Sato, M. 3D Ground Penetrating Radar to Detect Tree Roots and Estimate Root Biomass in the Field. Remote Sens. 2014, 6, 5754–5773. [Google Scholar] [CrossRef] [Green Version]
  28. Grasmueck, M.; Viggiano, D.A. Integration of ground-penetrating radar and laser position sensors for real-time 3-D data fusion. IEEE Trans. Geosci. Remote Sens. 2007, 45, 130–137. [Google Scholar] [CrossRef]
  29. Böniger, U.; Tronicke, J. On the potential of kinematic GPR surveying using a self-tracking total station: Evaluating system crosstalk and latency. IEEE Trans. Geosci. Remote Sens. 2010, 48, 3792–3798. [Google Scholar] [CrossRef]
  30. Rial, F.I.; Pereira, M.; Lorenzo, H.; Arias, P.; Novo, A. Resolution of GPR bowtie antennas: An experimental approach. J. Appl. Geophys. 2009, 67, 367–373. [Google Scholar] [CrossRef]
  31. Böniger, U.; Tronicke, J. Improving the interpretability of 3D GPR data using target-specific attributes: Application to tomb detection. J. Archaeol. Sci. 2010, 37, 672–679. [Google Scholar] [CrossRef]
  32. Bano, M. Modelling of GPR waves for lossy media obeying a complex power law of frequency for dielectric permittivity. Geophys. Prospect. 2004, 52, 11–26. [Google Scholar] [CrossRef] [Green Version]
  33. Hanafy, S.; Al Hagrey, S.A. Ground-penetrating radar tomography for soil-moisture heterogeneity. Geophysics 2005, 71, K9–K18. [Google Scholar] [CrossRef]
  34. Misra, P.; Enge, P. Special issue on global positioning system. Proc. IEEE. 1999, 87, 3–15. [Google Scholar] [CrossRef]
  35. Dow, J.M.; Neilan, R.E.; Rizos, C. The international GNSS service in a changing landscape of global navigation satellite systems. J. Geod. 2009, 83, 191–198. [Google Scholar] [CrossRef]
  36. Groves, P.D. Principles of GNSS, Inertial, and Multisensor Integrated Navigation Systems; Artech House: Norwood, MA, USA, 2013. [Google Scholar]
  37. Hofmann-Wellenhof, B.; Lichtenegger, H.; Wasle, E. GNSS–Global Navigation Satellite Systems: GPS, GLONASS, Galileo, and More; Springer Science & Business Media: Berlin, Germany, 2007. [Google Scholar]
  38. Jin, S.; Feng, G.P.; Gleason, S. Remote sensing using GNSS signals: Current status and future directions. Adv. Space Res. 2011, 47, 1645–1653. [Google Scholar] [CrossRef]
  39. Hegarty, C.J.; Chatre, E. Evolution of the global navigation satellitesystem (gnss). Proc. IEEE 2008, 96, 1902–1917. [Google Scholar] [CrossRef]
  40. Teunissen, P.J. Integer least-squares theory for the GNSS compass. J. Geod. 2010, 84, 433–447. [Google Scholar] [CrossRef] [Green Version]
  41. Montenbruck, O.; Hauschild, A.; Steigenberger, P. Differential code bias estimation using multi-GNSS observations and global ionosphere maps. Navig. J. Inst. Navig. 2014, 61, 191–201. [Google Scholar] [CrossRef]
  42. Hoque, M.M.; Jakowski, N. Higher order ionospheric effects in precise GNSS positioning. J. Geod. 2007, 81, 259–268. [Google Scholar] [CrossRef]
  43. Jin, R.; Jin, S.; Feng, G. M_DCB: Matlab code for estimating GNSS satellite and receiver differential code biases. GPS Solut. 2012, 16, 541–548. [Google Scholar] [CrossRef]
  44. Closas, P.; Fernández-Prades, C.; Fernández-Rubio, J.A. Maximum likelihood estimation of position in GNSS. IEEE Signal Process. Lett. 2007, 14, 359–362. [Google Scholar] [CrossRef] [Green Version]
  45. Wang, N.; Yuan, Y.; Li, Z.; Montenbruck, O.; Tan, B. Determination of differential code biases with multi-GNSS observations. J. Geod. 2016, 90, 209–228. [Google Scholar] [CrossRef]
  46. Baselga, S.; García-Asenjo, L. GNSS differential positioning by robust estimation. J. Surv. Eng. 2008, 134, 21–25. [Google Scholar] [CrossRef]
  47. Pan, S.G.; Gao, W.; Wang, S.L.; Meng, X.; Wang, Q. Analysis of ill posedness in double differential ambiguity resolution of BDS. Surv. Rev. 2014, 46, 411–416. [Google Scholar] [CrossRef] [Green Version]
  48. Rife, J. Collaborative vision-integrated pseudorange error removal: Team-estimated differential GNSS corrections with no stationary reference receiver. IEEE Trans. Intell. Transp. Syst. 2011, 13, 15–24. [Google Scholar] [CrossRef]
  49. Petovello, M.G.; O’Keefe, K.; Lachapelle, G.; Cannon, M.E. Consideration of time-correlated errors in a Kalman filter applicable to GNSS. J. Geod. 2009, 83, 51–56. [Google Scholar] [CrossRef]
  50. Won, J.H.; Pany, T.; Eissfeller, B. Characteristics of Kalman filters for GNSS signal tracking loop. IEEE Trans. Aerosp. Electron. Syst. 2012, 48, 3671–3681. [Google Scholar] [CrossRef]
  51. Brodin, G.; Daly, P. GNSS code and carrier tracking in the presence of multipath. Int. J. Satellite Commun. 1997, 15, 25–34. [Google Scholar] [CrossRef]
  52. Choy, S.; Bisnath, S.; Rizos, C. Uncovering common misconceptions in GNSS Precise Point Positioning and its future prospect. GPS Solut. 2017, 21, 13–22. [Google Scholar] [CrossRef]
  53. Yigit, C.O. Experimental assessment of post-processed kinematic Precise Point Positioning method for structural health monitoring. Geomat. Nat. Hazards Risk 2016, 7, 360–383. [Google Scholar] [CrossRef] [Green Version]
  54. Realini, E.; Caldera, S.; Pertusini, L.; Sampietro, D. Precise Gnss Positioning Using Smart Devices. Sensors 2017, 17, 2434. [Google Scholar] [CrossRef] [Green Version]
  55. Grasmueck, M.; Weger, R.; Horstmeyer, H. Full-resolution 3D GPR imaging. Geophysics 2005, 70, K12–K19. [Google Scholar] [CrossRef]
  56. McClymont, A.F.; Green, A.G.; Streich, R.; Horstmeyer, H.; Tronicke, J.; Nobes, D.C.; Pettinga, J.; Campbell, J.; Langridge, R. Visualization of active faults using geometric attributes of 3D GPR data: An example from the Alpine Fault Zone, New Zealand. Geophysics 2008, 73, B11–B23. [Google Scholar] [CrossRef] [Green Version]
  57. Novo, A.; Lorenzo, H.; Rial, F.I.; Solla, M. From pseudo-3D to full-resolution GPR imaging of a complex Roman site. Near Surf. Geophys. 2012, 10, 11–15. [Google Scholar] [CrossRef]
  58. Zhuge, X.; Yarovoy, A.G.; Savelyev, T.; Ligthart, L. Modified Kirchhoff migration for UWB MIMO array-based radar imaging. IEEE Trans. Geosci. Remote Sens. 2010, 48, 2692–2703. [Google Scholar] [CrossRef]
Remotesensing 12 01300 g001 550
Figure 1. (a) Surveying scenario of the Metasequoia trees using the proposed 3D GPR system. (b) The low-cost GlobalSat GPS receiver used in this investigation.
Figure 1. (a) Surveying scenario of the Metasequoia trees using the proposed 3D GPR system. (b) The low-cost GlobalSat GPS receiver used in this investigation.
Remotesensing 12 01300 g001
Remotesensing 12 01300 g002 550
Figure 2. Time-domain reflectometry (TDR) measurements of the soil moisture for the estimation of the subsurface velocity in the site.
Figure 2. Time-domain reflectometry (TDR) measurements of the soil moisture for the estimation of the subsurface velocity in the site.
Remotesensing 12 01300 g002
Remotesensing 12 01300 g003 550
Figure 3. Raw radargram of the tree root investigation using a 500 MHz RAMAC shielded antenna system.
Figure 3. Raw radargram of the tree root investigation using a 500 MHz RAMAC shielded antenna system.
Remotesensing 12 01300 g003
Remotesensing 12 01300 g004 550
Figure 4. (a) The location of the GEONET stations in Japan; (b) location of the investigation site and nearest GEONET Station.
Figure 4. (a) The location of the GEONET stations in Japan; (b) location of the investigation site and nearest GEONET Station.
Remotesensing 12 01300 g004
Remotesensing 12 01300 g005 550
Figure 5. (a) GPS receiver geometry with a reference station; (b) flowchart of the GPS receiver bias estimation and removal by use of the Kalman filter.
Figure 5. (a) GPS receiver geometry with a reference station; (b) flowchart of the GPS receiver bias estimation and removal by use of the Kalman filter.
Remotesensing 12 01300 g005
Remotesensing 12 01300 g006 550
Figure 6. Moving trajectories of antennas recorded by the GlobalSat GPS receiver. (a) The record by the GlobalSat GPS receiver; (b) post-processing results.
Figure 6. Moving trajectories of antennas recorded by the GlobalSat GPS receiver. (a) The record by the GlobalSat GPS receiver; (b) post-processing results.
Remotesensing 12 01300 g006
Remotesensing 12 01300 g007 550
Figure 7. Migrated vertical profiles along the survey direction. (a) Migrated profile at 0.8 m cross-survey direction; (b) migrated profile at 1.6 m cross-survey direction; (c) migrated profile at 2.4 m cross-survey direction; (d) migrated profile at 3.2 m cross-survey direction.
Figure 7. Migrated vertical profiles along the survey direction. (a) Migrated profile at 0.8 m cross-survey direction; (b) migrated profile at 1.6 m cross-survey direction; (c) migrated profile at 2.4 m cross-survey direction; (d) migrated profile at 3.2 m cross-survey direction.
Remotesensing 12 01300 g007
Remotesensing 12 01300 g008 550
Figure 8. Migrated horizontal slices at different depths. Red lines indicate the detected roots in the migrated data set. (a) Migrated slice at 16-cm depth; (b) migrated slice at 19-cm depth; (c) migrated slice at 27.5-cm depth.
Figure 8. Migrated horizontal slices at different depths. Red lines indicate the detected roots in the migrated data set. (a) Migrated slice at 16-cm depth; (b) migrated slice at 19-cm depth; (c) migrated slice at 27.5-cm depth.
Remotesensing 12 01300 g008
Remotesensing 12 01300 g009 550
Figure 9. Migrated horizontal slices at different depths. Red lines indicate the detected roots in the migrated data set. (a) Migrated slice at 42-cm depth; (b) migrated slice at 59-cm depth; (c) migrated slice at 79.5-cm depth.
Figure 9. Migrated horizontal slices at different depths. Red lines indicate the detected roots in the migrated data set. (a) Migrated slice at 42-cm depth; (b) migrated slice at 59-cm depth; (c) migrated slice at 79.5-cm depth.
Remotesensing 12 01300 g009
Remotesensing 12 01300 g010 550
Figure 10. The coarse root detection in a 3D migrated cubic: (a) B-scan profile which contains a tree root. (b) Root extension tracking with the −3 dB positions of a local peak (survey direction). (c) Root extension tracking with the −3 dB positions of a local peak (cross-survey direction).
Figure 10. The coarse root detection in a 3D migrated cubic: (a) B-scan profile which contains a tree root. (b) Root extension tracking with the −3 dB positions of a local peak (survey direction). (c) Root extension tracking with the −3 dB positions of a local peak (cross-survey direction).
Remotesensing 12 01300 g010
Remotesensing 12 01300 g011 550
Figure 11. The reconstructed 3D root system: (a) View from the starting point; (b) view from the ending point.
Figure 11. The reconstructed 3D root system: (a) View from the starting point; (b) view from the ending point.
Remotesensing 12 01300 g011
Table
Table 1. System parameters of the RAMAC 500 MHz Shielded Antenna system used for investigation purposes [27].
Table 1. System parameters of the RAMAC 500 MHz Shielded Antenna system used for investigation purposes [27].
ParametersRAMAC 500 MHz Shielded Antenna
Limit frequency lower than −10 dB: F L o w 138 MHz
Limit frequency higher than −10 dB: F H i g h 591 MHz
Bandwidth of −10 dB: B − 10 d B 453 MHz
Center frequency: F C e n t e r 364 MHz

Share and Cite

MDPI and ACS Style

Zou, L.; Wang, Y.; Giannakis, I.; Tosti, F.; Alani, A.M.; Sato, M. Mapping and Assessment of Tree Roots Using Ground Penetrating Radar with Low-Cost GPS. Remote Sens. 2020, 12, 1300. https://doi.org/10.3390/rs12081300

AMA Style

Zou L, Wang Y, Giannakis I, Tosti F, Alani AM, Sato M. Mapping and Assessment of Tree Roots Using Ground Penetrating Radar with Low-Cost GPS. Remote Sensing. 2020; 12(8):1300. https://doi.org/10.3390/rs12081300

Chicago/Turabian Style

Zou, Lilong, Yan Wang, Iraklis Giannakis, Fabio Tosti, Amir M. Alani, and Motoyuki Sato. 2020. "Mapping and Assessment of Tree Roots Using Ground Penetrating Radar with Low-Cost GPS" Remote Sensing 12, no. 8: 1300. https://doi.org/10.3390/rs12081300

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop