A method for calculating field crop green area index based on voxel size optimization
By employing 3D laser scanning and voxel size optimization algorithms, the problem of obtaining the green area index of wheat has been solved, enabling rapid and accurate calculation of the green area index, which is applicable to crop monitoring under different conditions.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- NANJING AGRICULTURAL UNIVERSITY
- Filing Date
- 2024-04-11
- Publication Date
- 2026-06-26
AI Technical Summary
Existing technologies struggle to quickly and accurately obtain the green area index of cereal crops such as wheat. Traditional methods are labor-intensive, time-consuming, and susceptible to external environmental factors. Therefore, there is a lack of effective methods for the application of ground-based lidar in this field.
Crop point clouds were acquired using a RIEGL VZ-1000 3D laser scanner. Through data registration, denoising, and voxelization, the green area index was calculated using a voxel size optimization algorithm. By combining the point cloud intensity information difference and the voxel size optimization algorithm, the extraction of green component point clouds and the calculation of the area index were realized.
It enables rapid and accurate extraction of wheat green area index with small error, and is applicable to different growth stages, varieties and nitrogen application levels, providing a non-destructive and rapid method for obtaining crop green area index in the field.
Smart Images

Figure CN118314195B_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of non-destructive monitoring of crop life information in precision agriculture, and relates to a method for calculating the green area index of crops in the field based on voxel size optimization, and particularly to a method for calculating the green area index of crops under field conditions using ground-based lidar point clouds. Background Technology
[0002] For cereal crops such as wheat, the green area index, which includes green organs such as leaves, stems, and ears, is an important indicator reflecting wheat growth status, yield, and quality. The green area index not only affects light distribution, light interception, and photosynthesis within the canopy, but is also closely related to crop growth rate, nitrogen allocation, and responses to (a)biotic stress. Therefore, timely and accurate monitoring of changes in the wheat green area index is essential in actual production.
[0003] Currently, most methods for obtaining the wheat green area index rely on traditional manual destructive sampling. This method is labor-intensive, time-consuming, and highly subjective, making it unsuitable for quickly and accurately obtaining wheat growth information and unsuitable for acquiring large-scale farmland growth information in precision agriculture. Furthermore, traditional optical remote sensing suffers from decreased accuracy in estimating the green area index later in crop growth due to the saturation effect of optical signals, and spectral information is easily affected by external environmental factors such as soil background and lighting conditions.
[0004] Ground-based lidar, as a phenotypic analysis tool for rapidly and accurately acquiring point clouds of targets over large areas and extracting three-dimensional structures, has been successfully applied to obtain crop canopy height and aboveground biomass. However, there is a lack of methods for calculating the green area index of wheat using ground-based lidar, especially for the rapid and accurate extraction of the green area index of field crops. Summary of the Invention
[0005] The purpose of this invention is to provide a method for calculating the green area index of crops in the field based on voxel size optimization. By using lidar to obtain the point cloud of crops growing in the field and using a voxel size optimization algorithm, the green area index of crops can be calculated based on ground-based lidar, which can overcome the drawbacks of traditional green area index acquisition methods.
[0006] To achieve the above objectives, the technical solution adopted by this invention is: a method for calculating the green area index of field crops based on voxel size optimization, which includes the following steps:
[0007] Step 1: Collect crop green area index samples and point clouds;
[0008] Step 2: Perform data registration and denoising on the crop point cloud obtained in Step 1, and then use the intensity information difference to obtain the green component point cloud;
[0009] Step 3: Voxelize the point cloud of crop green components obtained in Step 2 to obtain a three-dimensional voxel model;
[0010] Step 4: Apply the voxel size optimization algorithm to obtain the crop green area index.
[0011] Furthermore, in step 1, a RIEGL VZ-1000 3D laser scanner is used to acquire crop point clouds and simultaneously acquire crop green area index samples.
[0012] Furthermore, the method for acquiring the crop point cloud is as follows: an 8-site lidar test scheme is adopted, the scanning mode is 60 mode, that is, the angular resolution is 0.06°, and one point data is acquired for every 0.06° rotation;
[0013] The method for obtaining the crop green area index samples is as follows: A crop test area of 2 rows × 0.2m is selected in each plot. The leaves, stems, and panicles after heading are separated, retaining only the green and healthy organ components. A leaf area meter LI 3000 is used to scan the separated green leaves to obtain the area of all green leaves. The diameter of the stems and panicles is measured using calipers, and the length is measured using a ruler. The area of the stems and panicles is calculated as if they were cylinders. Finally, the leaf area, stem area, and panicle area are added together and divided by the corresponding sampling area to calculate the green area index for each plot.
[0014] Furthermore, in step 2, point cloud registration includes: using the ICP algorithm in the point cloud processing software RiSCAN PRO to register different coordinate systems to the same coordinate system; for two point cloud datasets A and B with a certain overlap, where A is the reference point cloud dataset and point cloud B is the coordinate system to be registered, they need to be rotated and translated to the coordinate system where point cloud A is located.
[0015] Furthermore, in step 2, the noise reduction process includes:
[0016] (1) Remove drift points: Use the set elevation threshold to remove noise points floating in the air. After filtering by elevation, abnormal point cloud data can be clearly seen and deleted.
[0017] (2) Remove miscellaneous points: Use the Deviation value to remove deviation points. The Deviation value is an indicator that measures the change between the transmitted pulse waveform and the echo pulse waveform. When the Deviation value is 0, it indicates that the waveform is not distorted and the shape of the received pulse completely coincides with the shape of the transmitted pulse. When the Deviation value exceeds the set value, it indicates that the waveform is severely distorted and these point cloud data need to be deleted.
[0018] Furthermore, after point cloud registration and denoising, the green component point cloud with higher moisture content absorbs laser energy at 1550nm, while the non-green component point cloud reflects laser energy. Based on the difference in point cloud intensity information, the green component point cloud of the crop is obtained.
[0019] Furthermore, in step 3, point cloud voxelization is performed as follows: based on the extracted green component point cloud boundaries, using the minimum value of point cloud X, Y, and Z (X... min Y min Z min Starting from a point (Δi), the region is divided into i*j*k finite voxels of size (Δi)*(Δj)*(Δk) with a step size of voxel size. The voxel coordinates and voxel values corresponding to the point cloud in the voxel coordinate system are determined. The voxel value is determined by judging whether the voxel contains the number of point clouds. If the number of point clouds in the voxel is greater than or equal to 1, it means that the laser beam has been intercepted and is marked as 1. If the number of point clouds in the voxel is equal to 0, it means that the laser beam has penetrated and is marked as 0.
[0020] Furthermore, in step 4, the voxel size optimization algorithm is as follows: The voxel size is determined based on the average distance of the target object's point cloud. For each query point, all neighboring points are determined within a pre-set distance, and all distances from the query point to these neighboring points are calculated. The minimum value among these distances is stored as d. min After querying all points, calculate and store d. min The average value is used as the voxel size; the crop green area index is estimated by calculating the probability that the laser beam penetrates the canopy and comes into contact with vegetation elements.
[0021] Furthermore, the probability formula for calculating the laser beam penetrating the canopy and contacting vegetation elements is as follows:
[0022]
[0023] Wherein, GAI is the green area index, n is the number of horizontal layers, ΔH is the height of the horizontal layer, and nl(k) and np(k) represent the laser beams intercepted and penetrated in the k-th height layer, respectively.
[0024] Furthermore, the feasibility and accuracy of the test algorithm were verified using sample data of the green area index in independent years, and the results were comprehensively evaluated using root mean square error (RMSE) and relative root mean square error (RRMSE).
[0025] The beneficial effects of this invention are as follows: The method for calculating the green area index of field crops provided by this invention acquires point clouds of field crops using lidar and applies a voxel size optimization algorithm to calculate the green area index of field crops, achieving rapid and accurate extraction of the green area index of field crops. This invention is novel and easy to operate, providing a theoretical basis and technical support for the rapid and non-destructive extraction of the green area index of field crops.
[0026] The field crop green area index calculated by the method of the present invention was compared with the field measured green area index. The verification results showed that the root mean square error (RMSE) was 0.51 and the relative root mean square error (RRMSE) was 16.48%, which demonstrated the accuracy of the method of the present invention.
[0027] Furthermore, through field trials with different growth stages, varieties, planting densities, and nitrogen application levels, this invention demonstrates that the green area index calculation method has excellent universality and accuracy. At planting densities of 25cm–40cm row spacing and nitrogen fertilizer levels of 0–300kg / ha, the method of this invention can significantly improve the estimation accuracy of the crop green area index throughout the entire growth period. Attached Figure Description
[0028] Figure 1 This is a flowchart of the method for calculating the green area index of wheat in the field according to the present invention.
[0029] Figure 2 It is a division of the point cloud of green and non-green components of wheat.
[0030] Figure 3 It is the dot clouding of wheat green components.
[0031] Figure 4 This is a result verification chart. Detailed Implementation
[0032] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments.
[0033] like Figure 1 As shown, this embodiment provides a method for calculating the green area index of field crops based on voxel size optimization, which mainly includes the following steps:
[0034] S01: Obtain the point cloud of the study area;
[0035] S02: Extraction of point cloud of green components from crops;
[0036] S03: Voxelization of green component point cloud data;
[0037] S04: Voxel size optimization algorithm;
[0038] S05: Verify the accuracy and universality of the algorithm using actual measurement data.
[0039] The specific details of the above steps are as follows:
[0040] S01: The RIEGL VZ-1000 3D laser scanner acquires point clouds of wheat in the field and simultaneously acquires measured data of the green area index of wheat in the field.
[0041] Two wheat varieties: Yangmai 15 (V1) and Yangmai 16 (V2). Two density levels: D1 row spacing is 25cm (2.4×10). 6 Seedlings / ha); D2 row spacing is 40cm (1.5×10 6 Seedlings / ha). Three nitrogen fertilizer levels were applied: 0 kg / ha (N0), 150 kg / ha (N1), and 300 kg / ha (N2). Nitrogen, phosphorus, and potassium fertilizers were urea, superphosphate, and potassium chloride, respectively. 50% of the nitrogen fertilizer was applied as basal fertilizer, 50% as topdressing at the jointing stage, and the remainder as basal fertilizer. Wheat was planted in rows. The experiment used a randomized block design, replicated three times, for a total of 36 plots. Each plot was 30 m². 2 (6m×5m), total area approximately 1080m² 2 .
[0042] Method for obtaining wheat point clouds: An 8-site LiDAR test scheme was adopted, with a scanning mode of 60°, i.e., an angular resolution of 0.06°. One point data was acquired for every 0.06° rotation.
[0043] Methods for obtaining field data on wheat green area index: For each plot, a wheat test area of 2 rows × 0.2m was selected. Wheat plant leaves, stems, and ears (after heading) were separated, retaining only the green and healthy organ components. A LI 3000 leaf area meter was used to scan the separated green leaves, obtaining the area of all green leaves. The diameter of the stems and ears was measured using calipers, and the length was measured using a ruler; the area of the stems and ears was calculated as if they were cylinders. Finally, the leaf area, stem area, and ear area were added together and divided by the corresponding sampling area to calculate the green area index for each plot.
[0044] Table 1. Specific details of the field trial
[0045]
[0046] S02: Perform data registration and noise reduction on the wheat point cloud obtained in step S01.
[0047] Point cloud registration: Each scanning station has its own independent coordinate system—the scanner coordinate system. In the point cloud processing software RiSCAN PRO, the ICP algorithm is used to register different coordinate systems to the same coordinate system (project coordinate system). The ICP algorithm is an optimal matching method based on least squares, repeatedly determining nearest neighbor pairs and calculating their coordinate transformation parameters until a given convergence accuracy is met, at which point the iteration ends. For two point cloud datasets A and B with some overlap, where A is the reference point cloud dataset and B is the coordinate system to be registered, rotation and translation transformation are required to bring them to the coordinate system of point cloud A.
[0048] Denoising: (1) Removing drift points: Noise points floating in the air are removed by setting an elevation threshold. After elevation filtering, abnormal point cloud data can be clearly seen and deleted. In this embodiment, the elevation threshold is set to 2m. (2) Removing mixed points: Deviation points are removed by setting a Deviation value. The Deviation value is an indicator that measures the change between the transmitted pulse waveform and the echo pulse waveform. When the Deviation value is 0, it means that the waveform is not distorted and the shape of the received pulse completely coincides with the shape of the transmitted pulse. When the Deviation value is large, it indicates that the waveform is severely distorted and these point cloud data need to be deleted. In this embodiment, the Deviation is set to 250.
[0049] After point cloud registration and denoising, based on the fact that green point clouds with higher moisture content absorb laser energy at 1550nm while non-green point clouds reflect laser energy, the difference in point cloud intensity information is utilized ( Figure 2 ), and obtained the point cloud of wheat green components.
[0050] S03: Based on the extracted green component point cloud boundary, using the minimum value of point cloud X, Y, Z (X min Y min Z min Starting from a point (Δi), and using the voxel size as the step size, the region is divided into i*j*k finite voxels of size (Δi)*(Δj)*(Δk). Figure 3 The process involves determining the voxel coordinates and voxel values of the point cloud in the voxel coordinate system. The voxel value is determined by checking the number of point clouds contained within the voxel. If the number of point clouds within the voxel is greater than or equal to 1, it indicates that the laser beam has been intercepted, marked as 1; if the number of point clouds within the voxel is 0, it indicates that the laser beam has penetrated, marked as 0.
[0051] S04: Voxel Size Optimization Algorithm
[0052] The voxel size is determined based on the average distance of the target object's point cloud. For each query point, all neighboring points are identified within a pre-set distance, and all distances from the query point to these neighboring points are calculated. The minimum value among these distances is stored as dmin. After querying all points, the average of the stored dmin values is calculated as the voxel size. The wheat green area index is estimated by calculating the probability that a laser beam penetrates the canopy and contacts vegetation elements.
[0053]
[0054] Wherein, GAI is the green area index, n is the number of horizontal layers, ΔH is the height of the horizontal layer, and nl(k) and np(k) represent the laser beams intercepted and penetrated in the k-th height layer, respectively.
[0055] S05: The accuracy and universality of the method in this embodiment are verified using wheat green area index data obtained from field tests, and the root mean square error (RMSE) and relative root mean square error (RRMSE) are used for comprehensive evaluation.
[0056] The calculation formulas for the evaluation indicators RMSE and RRMSE are as follows:
[0057]
[0058]
[0059] Where n is the number of samples, P i O is the model's predicted value. i These are experimental observation values. To observe the average value.
[0060] like Figure 4 The verification results showed that the root mean square error (RMSE) was 0.51 and the relative root mean square error (RRMSE) was 16.48%.
[0061] The foregoing has shown and described the basic principles, main features, and advantages of the present invention. Those skilled in the art should understand that the above embodiments do not limit the scope of protection of the present invention in any way, and all technical solutions obtained by equivalent substitution or other means fall within the scope of protection of the present invention.
[0062] All parts not covered in this invention are the same as or can be implemented using existing technologies.
Claims
1. A method for calculating the green area index of field crops based on voxel size optimization, characterized in that, Includes the following steps: Step 1: Collect crop green area index samples and point clouds; Step 2: Perform data registration and denoising on the crop point cloud obtained in Step 1. After point cloud registration and denoising, utilize the fact that the green component point cloud absorbs laser energy at 1550 nm, while the non-green component point cloud reflects laser energy. Based on the difference in point cloud intensity information, obtain the crop green component point cloud. Step 3: Voxelize the point cloud of crop green components obtained in Step 2 to obtain a three-dimensional voxel model; Step 4: Apply the voxel size optimization algorithm to obtain the crop green area index. The voxel size optimization algorithm is as follows: determine the voxel size based on the average distance of the target object's point cloud. For each query point, determine all neighboring points within a pre-set distance, calculate all distances from the query point to these neighboring points, and store the minimum value as d. min After querying all points, calculate and store d. min The average value is used as the voxel size; the crop green area index is estimated by calculating the probability that the laser beam penetrates the canopy and comes into contact with vegetation elements.
2. The method for calculating the green area index of field crops based on voxel size optimization according to claim 1, characterized in that, In step 1, a RIEGL VZ-1000 3D laser scanner is used to acquire crop point clouds and simultaneously acquire crop green area index samples.
3. A method for calculating the green area index of field crops based on voxel size optimization according to claim 1 or 2, characterized in that, The method for acquiring crop point clouds is as follows: an 8-site lidar test scheme is adopted, with a scanning mode of 60°, i.e., an angular resolution of 0.06°, and one point data is acquired for every 0.06° rotation. The method for obtaining the crop green area index samples is as follows: A crop test area of 2 rows × 0.2 m is selected in each plot. The leaves, stems, and panicles after heading are separated, retaining only the green and healthy organ components. A leaf area meter LI 3000 is used to scan the separated green leaves to obtain the area of all green leaves. The diameter of the stems and panicles is measured using calipers, and the length is measured using a ruler. The area of the stems and panicles is calculated as if they were cylinders. Finally, the leaf area, stem area, and panicle area are added together and divided by the corresponding sampling area to calculate the green area index for each plot.
4. The method for calculating the green area index of field crops based on voxel size optimization according to claim 1, characterized in that, In step 2, point cloud registration includes: using the ICP algorithm in the point cloud processing software RiSCAN PRO to register different coordinate systems to the same coordinate system; for two point cloud datasets A and B with overlapping parts, where A is the reference point cloud dataset and point cloud B is the coordinate system to be registered, it is necessary to rotate and translate them to the coordinate system where point cloud A is located.
5. The method for calculating the green area index of field crops based on voxel size optimization according to claim 4, characterized in that, In step 2, the noise reduction process includes: (1) Remove drift points: Use the set elevation threshold to remove noise points floating in the air. After elevation filtering, abnormal point cloud data can be clearly seen and deleted. (2) Remove mixed points: Use the Deviation value to remove deviation points. The Deviation value is an indicator that measures the change between the transmitted pulse waveform and the echo pulse waveform. When the Deviation value is 0, it indicates that the waveform is not distorted and the shape of the received pulse is completely consistent with the shape of the transmitted pulse. When the Deviation value exceeds the set value, it indicates that the waveform is severely distorted and these point cloud data need to be deleted.
6. The method for calculating the green area index of field crops based on voxel size optimization according to claim 1, characterized in that, In step 3, point cloud voxelization is performed as follows: based on the extracted green component point cloud boundaries, the minimum value of point cloud X, Y, and Z (X... min Y min Z min Starting from a voxel, the region is divided into i*j*k finite voxels of size (Δi)*(Δj)*(Δk) with a step size of voxel size. The voxel coordinates and voxel values corresponding to the point cloud in the voxel coordinate system are determined. The voxel value is determined by judging whether the voxel contains the number of point clouds. If the number of point clouds in the voxel is greater than or equal to 1, it means that the laser beam is intercepted and marked as 1; if the number of point clouds in the voxel is equal to 0, it means that the laser beam penetrates and marked as 0.
7. The method for calculating the green area index of field crops based on voxel size optimization according to claim 1, characterized in that, The formula for calculating the probability of a laser beam penetrating the canopy and contacting vegetation elements is as follows: SUBJECT = Wherein, GAI is the green area index, n is the number of horizontal layers, 𝛥𝐻 is the height of the horizontal layer, and 𝑛l(𝑘) and 𝑛𝑝(𝑘) represent the laser beams intercepted and penetrated in the k-th height layer, respectively.
8. The method for calculating the green area index of field crops based on voxel size optimization according to claim 1, characterized in that, The feasibility and accuracy of the test algorithm were verified using sample data of green area index from independent years. The results were comprehensively evaluated using root mean square error (RMSE) and relative root mean square error (RRMSE).