3.1 Data preprocessing
Different 3D sensing data need to be first processed into point clouds with different methods before sharing similar point processing methods. TLS data at different scanning locations were automatically registered to generate a point cloud using SCENE software (FARO Technology Inc, FL, USA). BLS was registered during data collection because the system used the SLAM algorithm (Su et al., 2021). GLS data registration was implemented according to the relative position of sensors and the point features using the commercial HortControl software (Phenospex Inc, Heerlen, The Netherlands). DAP images were used to reconstruct the 3D point cloud using the PiX4D mapper software (Pix4D, Lausanne, Switzerland). Once the 3D point clouds were generated, the following data processing processes were similar (Fig. 3).
The generated 3D point cloud data were further processed with a standard pipeline using the LiDAR360 software (Green Valley International Ltd., Beijing, China), including clipping, denoising, filtering, and normalization (Fig. 3). Clipping and denoising were manually implemented to ensure better accuracy, especially avoiding the loss of points in the sparse DAP and BLS point cloud. Filtering was first implemented using an integrated algorithm (i.e., improved progressive triangulated irregular network densification filtering algorithm), and the automatic results were carefully checked and revised to decrease process errors. Normalization was achieved by subtracting the height of each point from the height of its nearest ground point in the horizontal direction. Specifically, GLS data was filtered with a given height threshold of 0.28 m and normalized during data collection. The normalized 3D point clouds of TLS, BLS, GLS, and DAP were shown in Fig. 2b, d, f, g. Taking pre-processed data at the heading stage as an example, the point density of TLS data is the highest (929021.12 pts/m2), followed by GLS (697092.18 pts/m2), DAP (40051.30 pts/m2), and BLS (17761.30 pts/m2). Meanwhile, the final point resolution, denoted by the average adjacent point distance, from fine to coarse was GLS (1.07 mm), TLS (2.46 mm), DAP (12.73 mm), and BLS (15.02 mm) (Table 2).
Table 2
Key information about the data quality of the preprocessed point clouds (taking data at the heading stage as an example) and the roughly estimated platform cost and data cost.
Data sources
|
Point density (pts/m2)
|
Point resolution (mm)
|
Data volume, GB
|
Platform cost,$
|
Data Cost, h
|
Collection
|
Preprocessing
|
Total
|
TLS
|
929021.12
|
2.46
|
11.40
|
46010.00
|
5.00
|
79.20
|
84.20
|
BLS
|
17761.30
|
15.02
|
0.22
|
70515.00
|
0.30
|
2.30
|
2.60
|
GLS
|
697092.18
|
1.07
|
8.68
|
1567000.00
|
8.00
|
2.00
|
6.50
|
DAP
|
40051.30
|
12.73
|
0.48
|
1253.60
|
0.50
|
15.70
|
16.20
|
Plot extraction is the prerequisite for CH extraction of each plot. Because different sources of point clouds have their sensor coordinate systems, this study manually aligned these data into the same coordinate origin and north-south directions in LiDAR360 software. After that, 480 plots of different source data at each growth stage can be extracted using a shared plot bounding box map defined manually (Fig. 3).
3.2 Canopy height extraction
CH can be extracted from the normalized point cloud using different statistical metrics. In this study, Hmax, the maximum z value of all normalized points, was extracted. Meanwhile, difference height quantiles from 99% quantile height (i.e., H99) to 80% quantile height (i.e., H80) with an interval of 1% were also extracted (Jin et al., 2019). These different height representations are compared and the optimal one was selected for comparing different sensing technologies.
3.3 Cross-comparisons of canopy height estimates from field measurement and 3D sensing
The accuracies of the CH measured by different 3D sensing data were compared with the field measurement, and the cross-comparisons of different 3D sensing performances were also evaluated. Specifically, the comparisons between sensor data with field measurement include TLS vs.FM, BLS vs.FM, GLS vs.FM, DAP vs.FM, and the cross-comparisons include TLS vs. BLS, BLS vs. DAP, DAP vs. TLS, TLS vs. GLS, BLS vs. GLS, and DAP vs. GLS.
This study further evaluated the accuracy of different methods with respect to different field-measured CH groups, LAI groups, and GS groups, which are important indicators of canopy structure (Luo et al., 2015; Ma and Liang, 2022) and affect the accuracy of CH monitoring. Four CH groups were considered, including 0.3–0.6 m (CH1), 0.6–0.8 m (CH2), 0.8-1 m (CH3), and 1-1.4 m (CH4). Each height group contains 360, 918, 501, and 141 plots, respectively. Four LAI groups were separated at 0–2 m2/m2 (LAI1), 2–4 m2/m2 (LAI2), 4–6 m2/m2 (LAI3), and 6–8 m2/m2 (LAI4). Each group contains 874, 641, 340, and 65 plots, respectively. Four compared growth stages were jointing stages, heading stages, flowering stages, and maturity stages.
Specifically, considering the scanning range and height threshold setting in filtering, the effective maximum height of the GLS system is 0.82 m. Therefore, only the plots that have a maximum measured height lower than 0.82 m were selected for comparison with GLS (1365 plots) in this study. Because there are a few plots belonging to the CH3 group and no plots belonging to the CH4 group, we only evaluated the GLS accuracies of CH1 and CH2 (360 and 918 plots, respectively).
The accuracy between the two compared groups was evaluated by Pearson’s correlation coefficient (r), root mean square error (RMSE), relative RMSE (RMSE%), Bias, and relative Bias (Bias%).
\(r=\sqrt{1-\frac{\sum {\left({y}_{i}-{\widehat{y}}_{i}\right)}^{2}}{\sum {\left({y}_{i}-\overline{{y}_{i}}\right)}^{2}}}\)
(1)
\(RMSE=\sqrt{\frac{1}{n}\sum _{i=1}^{n}{({y}_{i}-{\widehat{y}}_{i})}^{2}}\)
(2)
\(RMSE\%=\left(\frac{RMSE}{\overline{{y}_{i}}}\right)\times 100\)
(3)
\(Bias={\sum }_{i=1}^{n}\left({y}_{i}-{\widehat{y}}_{i}\right)/n\)
(4)
\(Bias\%=\left(\frac{Bias}{\overline{{y}_{i}}}\right)\times 100\)
(5)
where i represents a sample index, n represents the number of samples, yi represents reference measurements (e.g., FM), \(\widehat{{y}_{i}}\) represents predicted CH from different 3D sensing datasets, and \(\stackrel{-}{{y}_{i}}\) is the mean of yi .
Moreover, the CHs of different data sources were compared in terms of broad-sense heritability (H2). Broad-sense heritability was defined as the proportion of heritability variance (Visscher et al., 2008). In this study, the interaction effect of different varieties and N treatments was considered, i.e., G by E.
$${H}^{2}=\frac{{\sigma }_{G}^{2}}{{\sigma }_{G}^{2}+\frac{{\sigma }_{GE}^{2}}{e}+\frac{{\sigma }_{\epsilon }^{2}}{re}}$$
6
where \({H}^{2}\) is broad-sense heritability, \({\sigma }_{G}^{2}\) is genetic variance,\({{\sigma }}_{\text{G}\text{E}}^{2}\) is the gene-environment interaction variance, \({{\sigma }}_{{\epsilon }}^{2}\) is the error variance, e is the number of N treatments, and r is the number of replicates per genotype.
3.4 Error source analysis
As we know, CHs measured by different methods will not be exactly the same. This study analyzed which data source the error comes from by referring to the method of Wang et al. (2019). First, we calculate the relative residual between the 3D sensing estimated CHs and FM (Eq. 7). Then, screening out the plots where the above calculated relative residuals greater than 20% as the suspicious cases (S) (Eq. 8). The intersections of STLS, SBLS, SGLS, and SDAP were defined as the errors due to FM (Error_FM) (Eq. 9). Based on Error_FM, the intersection of STLS, SBLS, SGLS, SDAP, and non-Error_FM was defined as the errors due to TLS (Error_TLS), BLS (Error_BLS), GLS (Error_GLS), and DAP (Error_DAP), respectively (Eq. 10–13). Notably, when regarding TLS or any other 3D sensing datasets as the errors, it is not mean the other three 3D sensing datasets do not contain outliers because the conditions for Error_FM are very strict.
The relative residual \({\varDelta }_{\left(a,field\right)}^{i}\), Sa, and Error_TLS, Error_BLS, Error_GLS, and Error_DAP were defined as below:
\({\varDelta }_{\left(a,field\right)}^{i}=\left|{H}_{a}^{i}-{H}_{filed}^{i}\right|/{H}_{field}^{i}\)
(7)
\({S}_{a}=\left\{{P}^{i}|{\varDelta }_{\left(a,field\right)}^{i}\ge 0.2\right\}\)
(8)
\(Error\_FM=\left\{{P}^{i}|{S}_{TLS}\cap {S}_{BLS}\cap {S}_{GLS}\cap {S}_{DAP}\right\}\)
(9)
\(Error\_TLS=\left\{{P}^{i}|{S}_{TLS}\cap (!Error\_field)\right\}\)
(10)
\(Error\_BLS=\left\{{P}^{i}|{S}_{BLS}\cap (!Error\_field)\right\}\)
(11)
\(Error\_GLS=\left\{{P}^{i}|{S}_{GLS}\cap (!Error\_field)\right\}\)
(12)
\(Error\_DAP=\left\{{P}^{i}|{S}_{DAP}\cap (!Error\_field)\right\}\)
(13)
where i represents a sample index, \({\varDelta }_{\left(a,field\right)}^{i}\) is the relative residual, \({H}_{a}^{i}\) and Sa represents predicted CH and the suspicious cases from 3D sensing datasets, where a can be TLS, BLS, GLS, and DAP, Error_FM, Error_TLS, Error_BLS, Error_GLS, and Error_GLS represent the errors from FM, TLS, BLS, GLS, and DAP, respectively.