Influence of flight altitude and control points in the georeferencing of images obtained by unmanned aerial vehicle

ABSTRACT This study aimed to explore the influence of flight altitude, density, and distribution of ground control points (GCPs) on the digital terrain model (DTM) in surveys conducted by unmanned aerial vehicles (UAVs). A total of 144 photogrammetric projects consisting of 399 aerial photos were carried out in a 2 ha area. These photogrammetric projects involved six GCP distributions (edge, center, diagonal, parallel, stratified, and random), six GCP densities, and four flight altitudes (30, 60, 90, and 120 m). The response surface methodology was used to find interference factors and total root-mean-square error (RMSEt) as well. The 60 m flight altitude presented was the most efficient. Central GCP distribution was observed to have low precision. Using stratified and random edge distributions, 10 GCPs are recommended to achieve geometric precision below 0.07 m at any flight height. However, for studies requiring up to 0.07 m precision, the best distribution was parallel with 4 GCPs at any altitude. Diagonal positioning of the GCPs showed RMSEt values below 0.11 m with 4 GCPs at any altitude. A good distribution of GCPs was found to be important, but the density of GCPs per image was more relevant when obtaining a lower RMSEt.


Introduction
Precision agriculture is a useful model in the management of natural resources and improvement of modern agriculture. (Far & Rezaei-Moghaddam, 2018;Orozco & Ramírez, 2016). Among precision agriculture techniques, aerial remote sensing presents new methods of research and work optimization, capturing terrestrial features using unmanned aerial vehicles (UAVs). Some applications of UAVs in agriculture were presented in studies of variables related to nitrogen in corn (Corti et al., 2019), the evaluation of water stress in agriculture (Gago et al., 2015), and precision agriculture (Mogili & Deepak, 2018).
Image collection using UAVs and their photogrammetric applications offers the possibility to observe agricultural fields from a different point of view. Therefore, it is possible to observe some field aspects that are relatively invisible when monitored from the ground (Candiago et al., 2015;Polo et al., 2015;Rodríguez-Fernández et al., 2017). In addition, UAVs offer other advantages, such as flexibility in collecting images, improved spatial resolution, and control over temporal resolution.
With the advent and popularization of UAVs for agricultural applications, photogrammetry has garnered interest and become one of the most modern technologies in crop management. However, it should be considered that the orthomosaic generated by aerial images presents geometric errors, which can be attenuated according to the terrain slope, image overlap, crop type, and flight altitude. Parameters such as image overlap and flight altitude are considered essential to optimize UAV flight missions. Variations in flight altitude are required in order to fly efficiently and faster.
Flight time is considered an important parameter to define flight plans, in some cases even compromising the study. One of the restrictions of aerial sensing using UAVs is the flight range (Traub, 2011). Aside from the work described in Ma et al. (2013), Deery et al. (2014) and Erdelj, Saif, Natalizio et al. (2019) explained that electric UAVs are unable to operate for long periods due to limited battery capacity. Commercial rotary wingtype UAVs typically achieve 25-30 min flights, thereby limiting continuous operation and large-scale coverage.
The general objective of photogrammetry is to represent characteristics of a surface with reliability in terms of precision and accuracy. (Daakir, Pierrot-Deseilligny, Bosser et al., 2017;Jalandoni et al., 2018). Photogrammetric reconstruction can be performed using software that applies structure from motion (SfM) algorithms (Izumida et al., 2017).
Among the photogrammetric techniques based on images, structure from motion (SfM) is one of the most used (Rahaman & Champion, 2019). With 2D images, epipolar geometry is estimated with resource matching algorithms, for example, the SIFT scaleinvariant feature transform (SIFT) algorithm (Ding. et al., 2018). This technique can be used to estimate external orientation on images and 3D object geometry reconstruction (Brandolini & Patrucco, 2019). Some approaches in the literature emphasize the capacity of SfM in the generation of digital elevation models (DEMs) (Castillo et al., 2012).
On the other hand, research indicated that there may be systematic deformations in UAV images (Rosnell & Honkavaara, 2012). In UAV applications, cameras generally capture images vertically and move parallel to the floor. Near-parallel imaging conditions and inaccurate self-calibration of unknown radial distortion can produce distorted reconstruction results (Li et al., 2016). These errors are known as a central domain (doming effect) and interfere with the geometric quality of the orthomosaic (Javernick et al., 2014). The doming effect is a fundamental problem of digital surface model (DSM) generation by SfM analysis associated with sets of almost parallel images and inaccurate correction of radial distortion of the lens (James & Robson, 2014).
There are two paths to achieving high precision in aerial surveys, i.e., aircraft equipped with real-time kinematic (RTK) systems or ground control points (GCPs) (Chiang et al., 2012;Tsai & Lin, 2017). GCPs can be obtained using a topographic or geodetic survey of points, for example, by using a total station or a pair of high-precision global navigation satellite system (GNSS) receivers. This method can support the georeferencing of data and geometric correction of images captured by UAVs, which are characterized as clearly visible reference targets in aerial images (Agüera-Vega, Carvajal-Ramirez, Martínez-Carricondo et al., 2018).
Studies such as those performed by Brunier et al. (2016) with real-time kinematic differential GPS (RTK-DGPS) and a total station generally prove to be of high precision and accuracy. However, this method requires a high spatial density of points to construct the digital elevation model (DEM); therefore, this type of work requires time, which increases the cost of the project. Ribeiro-Gomes et al. (2016) reported that the manual tasks needed to generate geometric products control the price of projects. Tasks such as walking in a field to collect points can slow a survey and require more than one operational professional, making the project more expensive.
The tracking system in a UAV consists of a GNSS receiver that provides its absolute location in the SIRGAS 2000 system. Because it is a navigation receiver, it produces a minimum error, in the order of meters (Zhang & Hsu, 2018). Errors in the range of meters are expected in GNSS receivers that collect only L1 data frequency. Differential GNSS (DGNSS) systems receive signals via two antennae of L1 and L2 frequency, thus collecting data from more satellites, increasing the number of triangulations, and improving precision and accuracy (Pervan et al., 2003).
Although control points bring geometric quality to the orthomosaic, some GCP investigations remain somewhat controversial, causing few studies to use this technique. In some cases, incorrect distribution and quantity of GCPs are used, i.e., the same precision could be achieved with a reduced number of points combined with adequate flight planning taking into consideration the photo scale, camera calibration, and flight project specification. Few studies in the literature correlate the influence of UAV flight altitude and the distribution and density of GCPs.
The time required for fieldwork is significantly optimized by the reduction of the number of GCPs (Eisenbeiss & Sauerbier, 2011). To increase the reliability of some variables such as flight efficiency and geometric precision in the orthomosaic it is essential to know the GCP distribution associated with the flight parameters. An efficient flight captures more images in less flight time maintaining acceptable geometric errors; this could decrease fieldwork time and accelerate project execution.
Therefore, this study aimed to evaluate the influence of the density of GCPs allocated in different distributions and submitted to different flight altitudes to assess the geometric precision, flight efficiency, and number of images in surveys conducted by UAVs.

Materials and methods
The methodology used in this work involves data collection steps via UAV and GNNS receivers, data processing, and interpretation of results. The processes of obtaining results are shown in flowchart form ( Figure 1) and detailed in the following sections.

Study area
The study was conducted in an experimental area consisting of coffee crops at the Federal University of Lavras (UFLA) (Figure 2), located in the municipality of Lavras, state of Minas Gerais, Brazil, covering an area of 2 ha with the geographical coordinates 21° 13ʹ33.23" south and 44°58ʹ17.63" west.

Acquisition of georeferenced data
The equipment used for data collection in the field is presented in Figure 3. In this area, 43 GCPs, predefined in a grid with 25 × 25 m intersections, were georeferenced with high precision and fixed in the field. The points were accurately obtained with an error of less than 0.03 m using a pair (base and rover) of GNSS antennas ( Figure 3A), Spectra Precision model SP60, operating in a real-time kinematic model (RTK).
To characterize and identify the GCPs, 0.3 × 0.3 m targets as shown in Figure 3C were placed at each tracked point obtained through the GNSS receivers. These GCPs were used to georeference the images obtained by the UAV.
EZSurv software and a digital platform of the Brazilian Institute of Geography and Statistics (IBGE) were used in order to process the collected points. This process consisted of transforming signals received between satellites and GNSS equipment into coordinates. To improve precision, the geographic coordinates (X, Y, and Z) from the RTK base were sent to the IBGE and adjusted by precise point positioning (PPP). Table 1 presents coordinates processed in the system of PPP. This positioning method applies an orbit and clock correction in the GNSS.  Its main benefit concerning differential positioning techniques is its ability to provide a position within a global reference frame anywhere in the world with a single GNSS receiver (Grinter & Roberts, 2011).
The coordinates obtained with the GNSS in rover mode were processed after processing the base coordinates, as shown in Figure 4. At this stage, the data between the GNSS base and GNSS rover were aligned in order to adjust the collected data precision. EZSurv   eliminated faulty incoming signals to achieve greater precision.

Acquisition of photogrammetric data
The aerial images were obtained with a DJI Phantom 4 Advanced aircraft as shown in Figure 3B, with the following characteristics: weight, 1,388 g; size, 350 mm; maximum speed, 72 km/h; maximum angle of inclination, 42°; maximum flight time, 30 min. The GPS/GLONASS positioning system was equipped, by default, with a 1 inch CMOS sensor to capture video (up to 4,096 × 2,160 p at 60 fps) and photos up to 20 megapixels. With the targets positioned at the tracking sites, the flight missions began. The planning was performed using the free software Drone Deploy installed on an Android 6.0 system. The flight plan was defined according to the following characteristics: speed, 3 m/s; front and side overlap, 60% × 80%, respectively; flight direction bearing, 50°; area, 2 hectares. The same mission was applied to the different altitudes of 30, 60, 90, and 120 meters.
Six different CGP distributions were defined: random, edge, center, diagonal, parallel, and stratified, as shown in Figure 5. The distributions were combined with variations in GCP density and different flight altitudes, forming a final combination that evaluated six distributions of 4, 5, 8, 10, 14, and 20 GCPs and flight altitudes of 30, 60, 90, and 120 m, totaling 399 aerial photos and 144 photogrammetric projects.
Processing of Photogrammetric Data Photogrammetric processing was performed using Agisoft PhotoScan software, version 1.4.3. According to Sona et al. (2014), this software, which is based on an SfM algorithm, is superior to others in terms of precision. It uses multiple camera views to increase photogrammetric data accuracy, not so different from aerial or terrestrial LiDAR. It can provide threedimensional points and produces a reliable data set to create dense point clouds. The input photographs can then be mosaicked and orthorectified to create the DEM by converting the point clouds into vector mesh or raster digital elevation models (DEMs) (Dietrich, 2016).
The methodology adopted by Flynn and Chapra (2014) and Rusnák et al. (2018) was used to generate the orthomosaic in six steps. In step 1, the alignment of the images was performed using the phototriangulation process and generation of a sparse Figure 5. GCPs combination where the x-axis represents the GCPs density and the y-axis represents the GCPs distribution. point cloud, which defined the coordinate system of the terrain. In step 2, the sparse point cloud generated in the previous step was densified for a more detailed representation of the mapped area and was also referenced the SIRGAS 2000 Zone 23S local coordinate system. In step 3, a model was built that accurately represented the three-dimensional mapped terrain. Thus, it was possible to represent the digital surface model (DSM), and, after filtering the point cloud of the soil, it was possible to visualize the digital terrain model (DTM). In step 4, the texture was applied to the model obtained in the previous step to improve the visual appearance and distinction between objects.
Step 5 consisted of the creation of the DEM. The generated products were two-dimensional raster format representations of the DSM and DTM. Lastly, the orthomosaic was generated in step 6.

Precision assessment
The process of geometric correction in images using GCPs consisted of a two-dimensional transformation, in which the coordinates collected in the pixel were replaced with the coordinates of the points obtained by receivers (GNSS) in the field (Tawfeik et al., 2016).
The orthomosaics for the altitudes of 30, 60, 90, and 120 m were georeferenced by PhotoScan 1.4.3. In this phase, the processed points obtained by the GNSS were exported in a txt file and later loaded into the PhotoScan program, thereby allocating GCPs at each georeferenced target position, so there was a combination of images with the same coordinate system. In this process, unfocused images or poorly positioned images in the photogrammetric alignment were relocated with greater precision, and differences are presented as root-mean-square error (RMSE).
RMSE is commonly used to express numerical accuracy results. It has the advantage of presenting values of errors in some dimensions of the analyzed variable (Hallak & Pereira Filho, 2011). It was considered an accuracy check to generate the RMSE data. This approach consists of validating the adjusted coordinates using independent georeferenced points. Each independent coordinate served as support to the studied GCP. The values used to calculate the total root-mean-square error (RMSEt) were the summation of the RMSE coordinate axis (X, Y, and Z), obtained by PhotoScan, and the calculation is performed using Equation 1.

RMSE t ¼
ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi where n is the number of GCPs; X Oi , Y Oi , and Z Oi are the X, Y, and Z coordinates that were respectively measured in the DEM; and X GNSSi , Y GNSSi , and Z GNSSi are the X, Y, and Z coordinates that were respectively measured with the GNSS in the field. The RMSEt value was used to compare all 144 photogrammetric projects obtained from 399 aerial photos in this study.

Statistical analysis
A response surface was obtained to represent the experimental data statistically using OriginPro 17 software. Response surface methodology (RSM) is one of the most widely used multivariate techniques to optimize processes (Nasri & Mozafari, 2018;Ronix et al., 2017). Based on the adaptation of a polynomial model to the experimental data, it is possible to predict the responses for all possible combinations of factors within a chosen experimental group (Bezerra et al., 2008), thereby determining a regression model to optimize an output variable, which is influenced by independent variables (Behera et al., 2018).
The values of the response variable (error) were plotted according to the sources of variation (altitude and density of points) on a three-dimensional graph of X (number of points in each distribution), Y (flight altitude), and Z (RMSEt). Then, a response surface was calculated to better represent each "position", as well as the surface equation and its respective correlation values.

Results and discussion
The results of the initial post-flight parameters based on the different flight altitudes are provided in Table 2. For altitudes of 30 to 120 meters, the flight time ranged from 2 to 12 min. However, between the altitudes of 60, 90, and 120 m, the variation was only 2 min.
The existence of a linear correlation between flight altitude, number of images, and spatial resolution is shown in Table 2, confirming that the higher the number of images collected, the better the spatial resolution results. This phenomenon, which was observed by Mesas-Carrascosa et al. (2016) revealed that spatial resolution was directly related to flight altitude and could be predefined to achieve greater detail in orthomosaic images.
As presented in Table 2, the 60 m flight exhibited the highest efficiency among the altitudes assessed, collecting 86 images in 4 min. When compared to the other flights, this altitude exhibited significantly interesting values regarding flight time, number of images, and spatial resolution. This superior performance may be related to the focal length of the sensor, the size of the analyzed area, and the flying height, thereby promoting the best use of complete images. Perroy et al. (2017) revealed the importance of flight altitude on the identification of species of trees in a forest compared to field sampling, concluding that as a flight reached higher altitudes, the possibility of correct identification decreased. In 2016, Quirós and Khot developed research relating flight altitude to the precision of counting plants in nurseries, highlighting the importance of flight altitude in the quality of the images and pointing out that above 40 m altitude, considerable errors were present in the plant counting results.
When observing the error presented by the different distributions studied (distribution of GCPs and flight altitude) (Figure 6), the low interference of flight altitude concerning georeferencing errors is evident. This analysis made it possible to identify the distribution that least met the proposed objectives. Figure 6 shows that the central distribution differed from the others, displaying errors of up to 0.12 m with the 30 m flight, which was above the mean of residual errors found in the other distributions.
It is also possible to see in Figure 6 the low performance of the center distribution. It is worth mentioning that the increase in flight altitude favored the decrease of RMSEt for this distribution. This situation occurred due to the reduction in the number of images, consequently increasing the number of georeferenced points in the orthomosaic. Furthermore, it was observed that the increase in flight altitude promotes a reduction in RMSEt. Son et al. (2019) studied the flight parameters and GCPs in the geometric quality of the orthomosaic obtained by a UAV and found better values of RMSEt in flight altitudes between 80 and 150 m.
Opposite results of RMSEt values were found at different altitudes without the use of GCPs. Rossi et al. (2017) observed that RMSEt errors were caused by factors such as flight altitude, lack of oblique images, low-cost camera, or higher relief. Given this, there are some contradictions for different sensors and types of terrain; however, in this study, good relationships were found between flight altitude and GCP distribution ( Figure 6).
Based on the results presented in Figure 7, the density of GCPs per hectare could be defined. Between 14 and 20 GCPs, a low variation in RMSEt values was verified. It was possible to observe that above 14  CGPs, the RMSEt values reduced by only 2 mm. The stability between 14 and 20 points in this study area was verified, and millimeter differences were observed in this interval. According to these observations, it was determined that no more than 14 GCPs were needed for any distribution studied in this work; therefore, it was possible to define a maximum of seven GCPs per hectare in fieldwork for flights above 30 meters. In 2017, Agüera-Vega, o showed that horizontal and vertical precision improved as the density of GCPs increased and also found a limit of GCPs where the RMESt values were shown to stabilize.
Based on the results presented in Table 3, the desired precision may be inferred for any distribution between altitudes of 30 to 120 m, with GCPs varying between 4 and 20 points. RMSEt values can be estimated before executing a flight mission through multiple nonlinear regression (Equation 2), where Z is the desired RMSEt. Thus, the values of the parameters shown in Table 3 added in Equation 2 are represented in Figure 8 where X is the quantity of GCPs and Y represents the flight height.
The correlation coefficient values (R2) range between 0.72 and 0.86, thereby demonstrating the reliability of the presented calculations. Figure 8 shows the interactions between the independent variables, flight altitude, and density of GCPs, where the low RMSEt values are represented by cold colors (blue) and the highest error values are represented by warm colors (red). The X-axis shows the density of points for any distribution stabilized with 14 GCPs at approximately 7 GCPs per hectare, thereby validating the information presented in Figure 7. Figure 8 shows that the results related to flight altitude for the diagonal and edge distributions behaved similarly. As shown in Figure 5, the diagonal distribution covered much of the border regions, thus providing similar characteristics to the graph presented in Figure 8.
As shown in the surface response graphs in Figure  8, the distribution of points in the central region of the area resulted in the highest RMSEt values, with the error reaching 20 cm in the 30 m flight altitude.
When compared to the 60, 90, and 120 m flight altitudes, this error tended to decrease because, at higher altitudes, the sensor collected fewer images (Table 2), which increased the density of points per image and improved adjustments via GCPs.
As seen in Figure 8, the best RMSEt results were expressed in the parallel distribution, possibly because this distribution model georeferenced images from the center to the edge distribution in a constant manner, thus covering a greater number of targets per captured image. The response results of the parallel distribution behaved in a way that improved the RMSEt by increasing the flight altitude. Figure 8 shows that the RMSEt values ranged from 3 to 7 cm, respectively, between 30 and 120 meters.
This research showed that flight altitude influenced RMSEt values, thereby counter-corroborating the results of Gómez-Candón, De Castro, & López-Granados (2014), who concluded in their research that flight altitude was an important parameter to consider when acquiring images using an UAV. However, they found no differences in the RMSEt georeferencing in the orthomosaics created by UAVs between 30 and 100 m high. Figure 8 shows that some distributions presented similarity between the graphs. Distributions such as diagonal, parallel, and stratified presented smaller RMSEt values, potentially because of the distribution pattern, i.e., the spacing between each GCP followed a uniform alignment for each direction, allowing the GCPs to be well distributed in the area. The influence of a poor distribution of GCPs was described in the research of Sanz-Ablanedo et al. (2018), who demonstrated that GCPs should be evenly distributed throughout the area of interest, ideally in an angular grid, because any GCP is minimized with maximum distance. These results indicated that, for a given density of GCPs, the precision obtained using optimal distribution would be twice as good as when GCPs were poorly distributed.
In the work of Martínez-Carricondo, Martínez-Carricondo et al. (2018), the authors suggested that to achieve more accurate values in orthomosaic georeferencing, GCPs should be distributed along the study area border in a stratified manner. However, Figure 8 shows a better distribution of the GCPs whereby parallel distribution was the most suitable for generating smaller error, i.e., below 7 cm and with a correlation coefficient of R2 = 0.83, potentially due to a more uniform distribution among the models considered.
The center, edge and random models exhibited higher RMSEt results, possibly due to the poor distribution of the GCPs and fewer points per image. Figure 8 demonstrates the low performance of the center distribution; at any flight altitude, this type of distribution georeferenced fewer images, as the number of points were clustered in the middle of the area.   (Table 2). Therefore, the greater points density per image, the better the geometric precision. For example, in the central distribution, which was considered the least precise distribution, as the altitude increased the number of GCPs increased as well, and the error decreased. This is because, in high flight altitudes for this distribution, it is possible to observe GCPs in the image center and also on the edge of the image. This was shown by the 120 m altitude flights.
As noted in Figure 9, the number of points per image tended to stabilize above 12 GCPs or six GCPs per hectare, or when the flight altitude exceeded 100 m, thereby corroborating the above discussion regarding the number of points per image. In the center distribution, as the number of points increased, the errors decreased, which caused the additional points to reach other images. Considering a spacing of 20 m at an altitude of 120 m, the center points may have georeferenced the edge images, which did not occur at low altitudes such as 30 m.
Mission planning prior to flight execution and field checkpoint surveys contribute to the objectivity of the work by reducing costs and operating time. In this research, the distribution of GCPs in the field and the previous adjustment of flight altitude were shown to contribute to a reduction in the number of GCPs and to achieving the previously defined RMSEt. Knowledge of the sensors attached to the aircraft enables more efficient missions regarding flight time.

Conclusion
The most efficient altitude was 60 m when considering the number of images and flight time. Based on the surface response analysis, the central distribution of GCPs was observed to have low precision. To achieve geometric precision below 0.07 m for stratified and random edge distributions, 10 GCPs at any flight height are recommended. However, for studies requiring up to 0.07 m precision, the best distribution was parallel with four GCPs at any altitude. The diagonal positioning of the GCPs showed RMSEt values below 0.11 m with four GCPs at any altitude. Despite good relationships between GCP density and distribution for the precision range, the number of image points had a greater influence on the georeferencing of the orthomosaic.