Introduction

Currently, photogrammetry is a technique used in obtaining digital elevation, and orthomosaic models of areas from aerial images (^{Mills & McLeod, 2013}). Although technology has potentially evolved to improve the results of this technique using tools based on terrestrial and remote sensors such as global positioning systems (GPS), total station, detection, and range of laser images (LiDAR), the synthetic aperture radar, and unmanned aerial vehicles (UAV), among others, the use of some of these continues to represent high costs in their use (^{Hastaoğlu et al., 2019}).

Among the instruments of high accessibility, low cost, and that allow obtaining high spatial resolution images in a simple, and flexible way, the UAV can be mentioned (^{Hemmelder et al., 2018}); this instruments allows safeguarding the integrity of operators in dangerous environments, and have few flight restrictions (^{ISPRS, 2004}; ^{Shi et al., 2011}; ^{Patterson et al., 2014}; ^{Xiang & Tian, 2011}).

Although UAVs have significantly increased their capacity to collect data in the field with higher density (^{Dering et al., 2019}), it has been pointed out that the use of additional tools improves the precision of the results (^{Jimenez et al., 2017}).

When only the camera positions are used as coordinates for the construction of the DEMs, these constructions tend to present errors in orientation of the block mainly, and open the way for there to be significant information biases. Therefore, it is recommended to add at least three or four GCPs to camera stations during processing to avoid bias because of elevations (^{Forlani et al., 2018}).

Thus, the use of Agisoft Meta Shape software and ground control points (GCP) added to the high-quality images contribute to improving the precision of the DEM obtained (^{Jimenez et al., 2017}).

The use of distributed ground control points (GCP), allows greater control in the plani-altimetry correction of the images, producing acceptable errors (^{Liew et al., 2012}; ^{Jimenez et al., 2017}; ^{Martinez et al., 2018}).

The software allows data processing and the generation of final models. However, these procedures make it possible to compare the methodologies with respect to the conventional methods. The bibliographic information is still scarce, and dispersed; the topography obtained by UAV requires case studies that contribute to validate the precision that can be achieved using this instrument, and the support of control points (^{Jimenez et al., 2017}).

Having said the above, this article describes a procedure to create a DEM and orthomosaic using a UAV and a GPS from high-resolution images contemplating eight georeferenced GCPs and corrected by the station closest to the study site.

Theoretical framework

Previous studies

Different authors have focused their work on the UAV and the GPS, due to the facilities and working conditions they provide, they also used these technologies for the development of their research such as ^{Niethammer et al. (2012)} and ^{Luccier et al. (2013)} used the UAVs to reconstruct landslides using aerial photography; ^{Siebert & Teizer (2014)}, they used UAV to inspect works and supervise the emergency and disaster service. ^{Barmpounakis et al. (2016)} used UAV to monitor traffic and acquire some parameters of speed, acceleration and trajectories. ^{Forlani et al. (2018)} used UAV with RTK (Real Time Kinematic or Real Time Satellite Navigation) in the same flight plan to generate Digital Surface Models (DSM) using three different block configurations.

Similarly, GPS has stood out as a new technology due to the ease of acquiring data in the field. In this way, ^{Llenera & Suarez (2015)} analyzed the GPS data to obtain reference and time systems using different positions of the satellites. They obtained points on the surface of the earth by means of computer triangulations, altitude, longitude, and latitude.

Now, referring to the study site that is proposed in this paper, some researchers have been carried out for some years where different damages caused by the existing faults in the Queretaro valley have been reported. These failures are mainly associated to the graben with shape rectangular bounded by the fault system with north-south direction and normal displacement. The western and eastern trace can be located in the San Bartolome fault and the Oriental fault (^{Carreón et al., 2005}). ^{Pacheco (2007)} carried out a gravimetric study, using six profiles on faulting, constructing gravimetric, and gradient anomaly maps in order to characterize the phenomenon on urban areas more precisely, quickly, and easy. On the other hand, ^{Civil Protection of Queretaro (2015)} created a Risk Atlas (Figure 1) in which it registers the vulnerability, risk of subsidence, and fracturing of the subsoil according to the location of each site, considering the proximity to the faults but this Atlas has an accuracy of 50 m in some areas and it does not detect the influence area of each failure nor the new ones. This Atlas is presented herein as Figure 1.

Figure 1 Sections with their respective GCP. Adapted from ^{Civil Protection (2015)}  

Materials and methods

Study area

The study area is located in Santiago de Queretaro, Queretaro, Mexico (Figure 2); it presents a rugged topography due to continuous movements, since there is a fault due to subsidence, also known by other authors as a stress crack, located in the Queretaro valley, which is shown in Figure 3 and Figure 4, crossing the neighborhood of "El Ensueño".

Figure 2 Macro and micro location of Santiago de Queretaro 

Figure 3 Fault system in the Central Zone and location in the “El Ensueño” (^{Protección Civil, 2015}) 

Figure 4 Micro-localization of the area of study “El Ensueño”, where the fault is shown by subsidence that crosses it 

This area presents vertical displacements of between 0.20 to 0.90 centimeters with NE-SW and NW-SE orientation. One of the most affected areas is the “El Ensueño” area, which is located inside the ellipsoid (Figure 5) WGS 84 (World Geodetic System, 1984) with the parameters of semi-major axis 6378137.00 and semi-minor axis 6356752.314, and is located in the zone 14 with UTM coordinates 2276546.01 N and 352689.16 E, where an area of 18,288 m² is considered (shown in Figure 3, inside a box).

Figure 5 Reference ellipsoid. Adapted from (^{Olaya, 2014}) 

Placements of ground control points

Eight GCPs were located at the study area in order to have greater control over the planimetry and altimetry correction of the images. GCPs were made up of steel plates welded to rods with an iron plate over it. This is shown in Figure 6. So, they were anchored to the ground and could be measured with a dual-frequency GPS (Figure 7).

Figure 6 Placement of GCP in the “El Ensueño 

Figure 7 Final result of the placement of the GCP 

As a result of the placement of the GCPs in the “El Ensueño” neighborhood, the following location sketch of the eight GCPs used is obtained, shown in Figure 8.

Figure 8 Location of the GCPs in the “El Ensueño” 

Furthermore, considering that the GCPs would not be visible during the UAV flight, the use of tarps with a multi-colored target was implemented to facilitate the identification of the GCPs.

Measurement of ground control points

The GCPs were measured on November 21 and 22, 2019, for which the Sokkia GRX2 GPS was used in static mode (Figures 9 and 10), in each of the eight points placed at the study site. Each point was measured for one hour in order to have better data quality and more reliable precision. So, the measurement error can be reduced.

Figure 9 GCP and GRX2 GPS used 

Figure 10 GCP measurement in the “El Ensueño” with the Sokkia GRX2 GPS 

Flight planning and image acquisition

The UAV flight plan was carried out on November 26, 2019 at 12:40 pm, in a time of 6 min and 26 s, with dimensions of 127 x 144 m and an overlap of images of 80 % and 70 %, at an altitude of 60 m. The flight plan carried out by the UAV at the study site is shown in Figure 11, all through the PIX4D mobile application.

Figure 11 UAV flight plan in the “El Ensueño” 

Images

For the processing of the images obtained by means of the UAV, the Agisoft MetaShape software was used where, firstly, the photographs obtained with the UAV were added, followed by the eight GCP obtained by means of the GPS data, and these were introduced as markers points, the coordinates obtained from the processing of the points with Topcon Tools were established, although it is worth mentioning that these points were attached manually, identifying in each photograph included in the processing, the GCP to have a minor error in the process; after that, the option of the cameras in each of the photographs was immediately unmarked, so that the coordinates for the model are those of the GCP and not those of the GPS of the UAV, because the latter has with many variations due to changes in climatic and environmental conditions, obstructions and sensors, among others.

Subsequently, the images were oriented with the coordinates obtained; as a result, the dense point cloud can be obtained, the model can be textured, the DEM and the orthomosaic can be acquired through the different integrated options that the software has.

Results

The results obtained from the processing of the data used in this article by means of GPS and UAV for obtaining DEM and orthomosaics are presented below.

The GCPs taken by the Sokkia GRX2 GPS in static mode, processed with Topcon Tools software and georeferenced to the INEGI RGNA, are shown in Figure 12 and Figure 13. The link of the active station of Queretaro and the GCPs are shown in UTM coordinates.

Figure 12 Processing of GCPs 1 to 4 with Topcon Tools software 

Figure 13 Processing of GCPs 5 to 8 with Topcon Tools software 

It should be noted that the processing was done separately, because the DEMO version of the software was used as it did not have the license.

As a result, the coordinates and elevations of the GCPs were obtained, such information is presented in Table 2.

On the other hand, the images that were obtained by the UAV were processed with the Agisoft MetaShape software, through which a dense point cloud anchored to the GCP obtained by the GPS was created for more precision of the model, finally acquiring the orthomosaic of the study site, which is shown in Figure 14.

Figure 14 Orthomosaic from the “El Ensueño” obtained by Agisoft Metashape Processing 

Likewise, the estimation of the GCP errors was obtained as shown in Figure 15. The report generated by Agisoft MetaShape where each color indicates the error in Z axis represented in the color palette on the right side, which are indicated in centimeters. The size and shape of the ellipse represents the error in X and Y, these errors are considered from the GCPs who are recognized as support points according to the software nomenclature. Therefore, when the images processed were obtained during the flight. These can be affected by the different light conditions and obstacles that the captured images may present.

Figure 15 Location of each GCP and error estimation using Agisoft MetaShape software 

Errors were automatically obtained by the software, using the formula for the calculation of the Mean Square Error (RMSE) by considering the square error in each of the dimensions with respect to its projection on the images that are shown in the previous in Figure 15. The GCPs presented in Figure 15 are shown below in Table 3.

In Table 3 each error for each GCP in the different dimensions are shown. The calculation of the errors were obtained by means of the differences between the source camera and the estimated coordinates of the GPS. By considering a good alignment of the images and the real coordinates obtained by the high-precision GPS. As a result of which a lower final value is obtained at 0.6 re-projection pixels for all points. It is considered that if the error obtained is greater than 1 pixel no correct alignment between the images and the coordinates is implied.

Finally, the digital elevation models of the study site are shown in Figure 16 and 17. A maximum height of 1,819.87 m and a minimum height of 1,795.59 m for the DSM were showed. A resolution of 0.0463611 m / pix, and a point density of 465,256 points per m², and for the DTM a maximum height of 1805.91 m and a minimum of 1789.45 m were obtained. The DEM is presented in Figure 16. It is better known as the Digital Surface Model (DSM) because it represents all the elements present on the study site or on the earth's surface such as: buildings, vegetation, automobiles, infrastructure, among others. Figure 17 is called the Digital Terrain Model (DTM). It shows a classification of the point cloud by terrain. Only those points that belong to the natural terrain, without considering all those points that belong to the infrastructure and vegetation are presented.

Figure 16 Reconstruction of the DSM 

Figure 17 Reconstruction of the DTM 

As can be seen in the Figure 17. The range of colors of blue varies in the central part of the DEM because in that area there is a subsidence transverse fault that completely crosses the “El Ensueño” neighborhood. So, it is possible to appreciate the unevenness by considering that it is approximately between 90 to 120 cm, which continuous displacements have been occurred causing different building’s damage during almost thirty years.

The parameters and the information mentioned above were obtained with the flight characteristics of the UAV shown in Table 4.

Conclusions

The UAVs have recently taken on an important role within engineering and research. Different authors have used the UAVs to carry out simple and sometimes more complicated tasks, some of them have been: topographic surveys by means of aerial photogrammetry, monitoring of large areas with landslides, emergency service in rescues from natural disasters, among other innovative uses. Due to the fact that their use is quite simple and they demand few resources for it use. It is being governed by the new (published past November 2019) established laws or regulations. Furthermore, they are considered an excellent alternative to acquire high resolution images and other data, in a short time in any environment.

On the other hand, GPS has also become an excellent tool for obtaining points at any physical site on earth, in order to accurately obtain geographic coordinates.

According to the aforementioned. It is exposed in the different sections of this article, the use of both technologies is complementary, increases precision and reduces data collection time in the field. However, the use of the UAV independently has some disadvantages when used to carry out photogrammetric surveys. Because the integration of GPS and UAVs have multiple variations due to the geometry of the satellite mainly due to the different weather conditions and the quality of the equipment.

As can be seen in Table 2, the elevations of the GCP were obtained quickly and easily, repeating this step for a certain time, a monitoring of the complete area could be carried out in order to observe and understand its behavior during related phenomena with displacements. On the other hand, in Table 3 the different errors can be observed for each GCP by considering the pixels of each one of the photographs use. As an example GCP called as BN-1 presents an error in centimeters at X = 4.44092, Y = 0.50457 and Z = 2.74891.

It must be mentioned that these errors are automatically calculated by the software using the formula of the root mean square error (RMSE). This depends on how a pixel is altered in the processing, or in the capture of the images. The advantage of this process through Agisoft MetaShape is that it shows the total errors in the X, Y and Z dimensions. It can be seen that these errors depends on the number of images to which they are georeferenced during post-processing. It can be seen that GPS and UAV are two potential elements for obtaining DEM and topography in a short period of time, and with low use of technical personnel to carry out measurements and operation of the flight, in addition to the low requirements that are needed for the processing of all the data.

The most important impact is that the Digital Terrain Model (DTM) shows a spectacular difference between one side and another of the failure zone, which can indicates that this is a perfect tool to distinguish the particular affected areas with to use in an Atlas risk for urban areas with construction over failures due to land subsidence phenomena were detected or not.