Coastal Mapping Using DJI Phantom 4 RTK in Post-Processing Kinematic Mode

: Topographic and geomorphological surveys of coastal areas usually require the aerial mapping of long and narrow sections of littoral. The georeferencing of photogrammetric models is generally based on the signalization and survey of Ground Control Points (GCPs), which are very time-consuming tasks. Direct georeferencing with high camera location accuracy due to on-board multi-frequency GNSS receivers can limit the need for GCPs. Recently, DJI has made available the Phantom 4 Real-Time Kinematic (RTK) (DJI-P4RTK), which combines the versatility and the ease of use of previous DJI Phantom models with the advantages of a multi-frequency on-board GNSS receiver. In this paper, we investigated the accuracy of both photogrammetric models and Digital Terrain Models (DTMs) generated in Agisoft Metashape from two different image datasets (nadiral and oblique) acquired by a DJI-P4RTK. Camera locations were computed with the Post-Processing Kinematic (PPK) of the Receiver Independent Exchange Format (RINEX) ﬁle recorded by the aircraft during ﬂight missions. A Continuously Operating Reference Station (CORS) located at a 15 km distance from the site was used for this task. The results highlighted that the oblique dataset produced very similar results, with GCPs (3D RMSE = 0.025 m) and without (3D RMSE = 0.028 m), while the nadiral dataset was affected more by the position and number of the GCPs (3D RMSE from 0.034 to 0.075 m). The introduction of a few oblique images into the nadiral dataset without any GCP improved the vertical accuracy of the model (Up RMSE from 0.052 to 0.025 m) and can represent a solution to speed up the image acquisition of nadiral datasets for PPK with the DJI-P4RTK and no GCPs. Moreover, the results of this research are compared to those obtained in RTK mode for the same datasets. The novelty of this research is the combination of a multitude of aspects regarding the DJI Phantom 4 RTK aircraft and the subsequent data processing strategies for assessing the quality of photogrammetric models, DTMs, and cross-section proﬁles.


Introduction
Nowadays, the need for mapping long coastal sections up to the survey of the littoral of entire regions (such as Emilia Romagna in Italy) represents the starting point for coastal management and for the planning of restoration interventions to mitigate the impact of storm events by regional environmental protection agencies.
Therefore, the production of Digital Terrain Models (DTMs) is a primary tool for the monitoring of coastal erosion, as well as for developing detailed and high-resolution numerical simulation models able to predict the impact of storms, sea level rise, and flooding risk [1]. Dune systems are characterized by high dynamics in which the morphology changes rapidly and where high accuracy is required [2] to describe the evolution of small dunes [3].
Among the multitude of surveying methodologies by which it is possible to adopt in coastal environments [4], the most important are certainly (i) the direct survey with GNSS receivers, (ii) Terrestrial Laser Scanning (TLS) [5], (iii) aerial survey by planes and drones with Light Detection And Ranging (LiDAR) [6], and (iv) the use of aerial photogrammetry [7]. However, direct methods are very time-consuming, and this makes them hard to apply to the mapping of large coastal section extents. Conversely, TLS is able to provide good spatial resolution, even if the duration of data acquisition is still high and the post-processing of scans generally requires an additional and considerable amount of time, while LiDAR by plane represents the most used technique for the mapping of Northern Adriatic Sea (Italy) in the last several years [8], because it ensures a reasonably quick data acquisition. In spite of this, the spatial resolution of final terrain models is usually no better than 1 m. Thus, the DTMs cannot be used for the high-resolution monitoring of limited extents of coastal sections up to a centimeter-level accuracy, nor can they be used for the planning of mitigation and restoration interventions on the beach. Recently, LiDAR systems on Unmanned Aerial Vehicles (UAVs) are becoming more accessible to users [9], but the cost of such platforms is still high in most cases. In addition, UAV LiDAR generally requires a considerable amount of time for the post-processing of point clouds.
The use of photogrammetric techniques allows the generation of dense point clouds [10] whose density is comparable to those acquired by TLS and LiDAR [11,12] and in which the Ground Sample Distance (GSD), for the same camera, depends on the flight altitude only. Traditional photogrammetry by plane generally requires the use of manned aerial vehicles and aerial cameras that are very expensive, and this makes it competitive only when a large extent of littoral has to be mapped. On the contrary, UAV photogrammetry is based on drones [13][14][15] at an affordable cost [16,17]. On one hand, the use of fixed-wing UAVs still enables the mapping of a substantially large beach area, but with operational limits in terms of the maximum distance allowed between the remote pilot and the aircraft depending on local regulations; on the other hand, multi-copters are extremely versatile and allow one to easily acquire oblique images due to the presence of a gimbal.
Substantially, surveying time using UAVs is related to the performance of flight missions initiated to acquire image datasets and for the deployment, survey, and recovery of Ground Control Points (GCPs). In fact, GCPs generally consist in targets that are useful (and usually needed) for an accurate georeferencing of the models, although operations related to GCPs during fieldwork are usually time-consuming. The geotagging of aerial images with standard GNSS data saves time during the identification of targets on the images, since the user can generate an approximate georeferenced model prior to any further operation related to GCPs. However, standard GNSS does not provide a sufficient constraint for the reconstruction of photogrammetric models with direct georeferencing due to the poor accuracy of exterior orientation parameters. In particular, camera locations have an uncertainty up to a certain amount of meters in the East (E), North (N), and Up (U) directions.
Whenever an aircraft is equipped with a multi-frequency GNSS receiver, the path followed by the drone during a flight is determinable a posteriori with a centimeter-level accuracy, and it is not necessary to use a large set of GCPs for model georeferencing, with a considerable reduction in the timing of survey. Thus, such drones are particularly interesting in the field of coastal mapping. Moreover, the acquisition of images can generate not only high-resolution DTMs but also orthomosaics at a spatial resolution of a few centimeters. The recent availability of UAVs with an on-board multi-frequency and multi-constellation GNSS receiver [18] makes direct georeferencing accessible to almost any kind of user at decreasing cost [19], and such aircrafts can record their own flight path through GNSS observations and geotag images with centimeter-level accurate coordinates. Strategies for using this technology are the same as for traditional geodetic GNSS receivers. Thus, Real-Time Kinematic (RTK) with a base and rover, Network RTK (NRTK) [20], and Post-Processing Kinematic (PPK) [21] of GNSS observations stored in a Receiver Independent Exchange Format (RINEX).
In coastal regions with low cellular network signals, receiving NRTK corrections may be particularly difficult. Conversely, the use of a base and rover requires the simultaneous use of two receivers with a considerable cost increase and requires one to put the base receiver on a benchmark, whose coordinates should be known in order to ensure the most accurate georeferencing of photogrammetric models. This issue can be overcome by performing a benchmark survey through NRTK if possible or, alternatively, with static GNSS. The PPK approach, instead, can be used whenever the aircraft can record a RINEX file with GNSS observations during a flight mission and when GNSS data acquired by local receivers or by a Continuously Operating Reference Station (CORS) are simultaneously available. In particular, the use of CORS data makes it possible to perform automatic flight missions without any need to use a base receiver or to connect to a network service for receiving NRTK corrections.
In this article, we investigated the PPK approach for the reconstruction of photogrammetric models and DTMs of a coastal section located in the Northern Adriatic Sea (Italy) through the use of RINEX files and the GNSS data of a CORS that is approximately 15 km away from our site ( Figure 1) using a recent multi-copter model (DJI Phantom 4 RTK). The aim of this research is to assess the quality of photogrammetric models and DTMs obtained through the PPK mode, whose processing strategies have already been proven to produce good results for the same datasets with RTK camera locations [22]. An evaluation of the real accuracy, which is possible to reach with such a multi-copter through post-processing with the GNSS data of a CORS at a distance of 15 km away from the site, is useful for investigating the mapping potential of the DJI Phantom 4 RTK of a littoral area where fixed-wing drones are usually used. In addition, after the analysis of the results of the PPK strategy, we also compared them with those obtained through RTK for the same aerial datasets. The novelty of this research is therefore the combination of all of these aspects regarding the DJI Phantom 4 RTK aircraft and the assessment of the results of different processing strategies, especially for what concerns the quality of photogrammetric models and DTMs, with a quick approach in which no geodetic GNSS receiver is ever used in the field and only images are captured, together with a RINEX file to be post-processed with CORS data.

Materials and Methods
The site mapped in this work is located along the coast of the Northern Adriatic Sea (Italy) near the city of Ravenna. The littoral consists in a coastal section with approximate dimensions of about 2 km in the North-South direction and a width of up to 130 m (see Figure 1). This coastal section has been selected for our research because it is very representative of the littoral of Emilia Romagna, since its narrow shape is very common in this region, with the only exception being beaches located near harbors. In addition, the possibility of having a mean CORS-drone distance of about 15 km was an added value of this study site for investigating the PPK approach.
The aerial surveys were performed using a recent model of a RTK multi-copter: the DJI Phantom 4 RTK (DJI-P4RTK), a compact and lightweight UAV with a 20 Mpix camera that has been available on the market since slightly more than one year ago. This model of aircraft for this research was used because it continues the great impact on users of previous Phantom DJI multi-copters, and it is also the most versatile and easy to use lightweight RTK DJI multi-copter.
The addition of a GNSS antenna+receiver module to enable the DJI-P4RTK to kinematic GNSS [23] with centimeter-level accuracy has therefore made this aircraft of great interest in the field of coastal mapping with direct georeferencing. In particular, the GNSS receiver, with which the DJI-P4RTK is equipped, is a multi-constellation multi-frequency one and is able to receive GPS, GLONASS, Galileo, and Beidou signals. The sampling rate of raw GNSS data is 5 Hz. The flight path is continuously recorded during each automatic flight mission, and the GNSS observations are stored into a RINEX file v.3.03. Whenever an image is captured by the camera of the DJI-P4RTK on a mission, a GPS timestamp (in terms of GPS Week and Time of Week) is also recorded.
Although all the procedures adopted in this work are presented in detail in the following paragraphs, a summary of the overall methodology and data processing is illustrated in Figure 2 for clarity.

In-Situ Operations
During flight operations, two different datasets of aerial images were acquired, nadiral and oblique. Moreover, in order to ensure a GSD of about 2 cm, all the flights were carried out at an altitude of 80 m. The flights were performed flying at 6 m/s, side and forward overlap were set to 70% and 80%, respectively, and the GS RTK app integrated into the remote controller of the DJI-P4RTK was used for flight planning.
The nadiral mapping was planned through an automatic single-grid mission, covering the entire site region and acquiring a total of 723 images. The flight plan for the acquisition of the oblique dataset has been limited to a smaller region in the center of the overall site. In this case, the flight plan consisted in a double-grid mission (camera was set with a pitch of 30°), and 612 images were acquired.
During fieldwork, about 40 points have been also signalized on the ground by a cross (Figure 3a). The cross size was designed for an easy detection on each acquired image. The points were deployed in order to be well-distributed along the entire area of the site and were further considered as GCPs for the georeferencing of the models or as Check Points (CPs) to assess both the horizontal and vertical accuracy of the models. The survey of GCPs/CPs has been performed through NRTK with 30 s stop-and-go. The network service broadcasting NRTK corrections used in this work is NetGEO (Topcon Positioning Italy).
Finally, more than 100 natural points were also surveyed on the beach, i.e., Validation Points (VPs), to evaluate and assess the vertical accuracy of DTMs (Figure 3b). Similarly to GCPs and CPs, all the VPs were surveyed through NRTK with 30 s stop-and-go, even if the survey was performed over the entire day and was not limited to a short time frame (i.e., just before the flights), as it was for the survey of GCPs and CPs.

Post-Processing Kinematic of RINEX Files
In this work, the RINEX file of each automatic flight mission has been post-processed using the GNSS observations of a CORS. The nearest NetGEO CORS was identified as the station whose ID is RAVE. This station is located in the city of Ravenna, and it is about 15 km away from our site. Furthermore, RAVE is currently able to receive only GPS and GLONASS signals with a sampling rate of 1 Hz. Therefore, no Galileo and Beidou data recorded by the DJI-P4RTK were used during PPK. The post-processing of flight paths has been carried out with Topcon Tools v.8.2.3 (a standard software for processing static and kinematic GNSS data) introducing RAVE coordinates in the official Italian reference system, corresponding to the European Terrestrial Reference System ETRS89 in its ETRF2000 Due to the GPS timestamp of each image, it has been possible to compute the APC position at the time of each image acquisition with an interpolation of the PPK flight path solutions along the E, N, and U directions (Figure 4b). Since the speed of the aircraft was 6 m/s, the timestamp should be measured with an accuracy better than 1.5 ms to ensure a centimeter-level accuracy of the final camera location due to the interpolation process itself. Timestamps are recorded by the DJI-P4RTK with a resolution of 1 µs.
Moreover, the offset between the APC and the camera center is known since it is reported among the technical specifications of DJI-P4RTK [24]. Referring to the axes shown in Figure 5 (fixed to the body of the aircraft), the x, y, z offsets are 36, 0, and 192 mm, respectively.
Together with the GPS time, offsets in the E, N, and U directions are also recorded during each image acquisition and stored within a text file. Using the E, N, and U offsets, it is possible to calculate the position of the camera center for every image captured during a flight mission simply by adjusting the corresponding APC coordinates (Figure 4c) [25,26]. Offsets along the E, N, and U directions are computed by the DJI-P4RTK autonomously during flight operations using the real-time yaw, pitch, and roll measured by the Inertial Measurement Unit (IMU).
A MATLAB script was then structured in order to give both coordinates and accuracies in the E, N, and U directions for each acquired image. The flight paths and the final camera locations of both datasets are shown in Figure 6. Minimum and maximum accuracy values, as they resulted after post-processing performed through Topcon Tools, are also reported in Table 1, together with the average.    Ellipsoidal heights of all the points (PPK camera locations, GCPs, CPs, and VPs) have been transformed into orthometric elevations by applying the ItalGEO 05 geoid separation model.

Generation of Photogrammetric Models and DTMs
After the determination of all the camera location coordinates through PPK, interpolation, and offset correction due to the APC position, photogrammetric models were generated for both nadiral and oblique datasets. For this task, we used Agisoft Metashape v.1.5.
In order to testing the potential of direct georeferencing with DJI-P4RTK with as few GCPs as possible, projects were created varying the number and the position of GCPs [27][28][29]. In particular, as already pointed out by Taddia et al. using RTK camera locations for the same datasets [22] as well as Forlani et al. [30] with a fixed wing aircraft, the nadiral dataset necessarily requires the use of at least one GCP for the model to obtain good accuracy. The lack of any GCP, indeed, may cause a wrong estimation of the focal length if only nadiral images are used, generating a model with a significant vertical offset of up to meters in some situations. Consequently, for the nadiral dataset, the case without the use of any GCP was not investigated here, but we began our analysis from the case with the introduction of a single GCP, either lateral or central with respect to the bounds of the flight path (Figure 7a). The case with the use of a set of well-distributed GCPs was also investigated (Figure 7b).
For the oblique dataset, we generated photogrammetric models for (i) the case without any GCP and (ii) using a set of GCPs. Furthermore, in this work, we constrained the camera locations to have an a priori accuracy equal to the standard deviation computed as shown in Section 2.2, while the other exterior orientation parameters (yaw, pitch, and roll) were left free to be adjusted by Agisoft Metashape (a default angular accuracy of 10°was set). Interior orientation parameters were always estimated through the self-calibration in this work.
Due to the large amount of images for both datasets, each of those with a resolution of 20 Mpix, we decided to perform the alignment with the parameter accuracy set to medium in Agisoft Metashape, and the following parameters were optimized: f , c x , c y , k 1−3 , p 1−2 , b 1 . The accuracy assessment of the alignment results was then made through a set of CPs. In particular, here we assumed as CPs all the crosses surveyed as described in Section 2.1 that were not used as GCP in a project.
For each alignment, we computed the Root Mean Square Error (RMSE) in the E, N, and U directions and in 3D as follows: where n represents the number of GCPs/CPs considered. After the assessment of the quality of the photogrammetric models on the basis of the CP residuals, computed as the difference between the position estimated through the model and the coordinates surveyed with GNSS NRTK, we proceeded with the densification of the point cloud (the densification level was set to medium with an aggressive filtering). Finally, we generated a DTM and evaluated the accuracy of it with a further quality assessment. For this task, we used the set of VPs and computed vertical differences between the DTM and the GNSS survey (∆H = H DTM − H GNSS ). We generated frequency histograms to evaluate the mean and standard deviation of the distribution of those differences.
After the processing of nadiral and oblique datasets, due to both the speed in acquiring images with single-grid nadiral flight missions for the mapping of large littoral areas as well as the need for introducing oblique images to estimate the focal length accurately, we created a "hybrid" data processing strategy combining the nadiral dataset with a short sequence of oblique images extracted from the second dataset. We selected 31 images. The aim was to increase the vertical accuracy of the DTM without any GCP. It is worth noting that the possibility of acquiring some oblique images during nadiral flight missions is also an option (called altitude optimization) of the GS RTK app integrated into the remote controller of the DJI-P4RTK. With this option enabled, the aircraft can capture a short sequence of oblique images from the position where the automatic mission ends and from the center of the area of the flight plan. Although we used this app for flight planning, we did not enable this option at the time of the aerial survey. However, this further "hybrid" data processing is also useful in simulating the results of overall image acquisition.

Nadiral Dataset
The results of the alignment of the nadiral image dataset confirmed the necessity of a set of well-distributed GCPs to guarantee the creation of photogrammetric models with a final centimeter-level accuracy. The use of a single GCP, either lateral or central with respect to the area covered by the aerial missions, can lead to vertical accuracies with differences of up to about 10 cm on CPs (Table 2). Apart from a vertical offset, the models in which a single GCP has been introduced show an RMSE that is comparable in the E, N, and U directions. Up RMSEs amounted to about 5 cm, while 3D RMSEs were ≈ 7.5 cm.
However, it is only by using a set of GCPs that the residuals decrease drastically, despite the use of PPK camera locations that are known with high accuracy. In this case, the maximum vertical CP residual is slightly less than 4 cm, with a 3D RMSE of 3.4 cm. Minimum, maximum, and RMSE values for CP residuals (along with GCP residuals for the case with a set of GCPs) for the nadiral dataset are reported in Table 2. Table 2. CP residuals for the nadiral dataset (GCP residuals in brackets). The subsequent validation of DTMs has been performed using 112 VPs for each photogrammetric processing strategy discussed above. In spite of a distribution of the elevation differences computed between the GNSS survey and the DTMs characterized by a standard deviation always around 4 cm, the mean value highlights a bias that is particularly evident when a GCP is assumed at the edge of the surveyed region (see Figure 8). Conversely, the cases in which the GCP was assumed in a central position or when a set of GCPs was adopted are located around a null average, deviating from it by about 2 cm.

Oblique Dataset
The results of the alignment of the oblique image dataset, evaluated in terms of residuals on the CPs, show that the use of a double-grid flight plan and a tilt of the camera makes it possible to reach a centimeter-level accuracy in the E, N, and U directions, even without the introduction of any GCP into the photogrammetric model.
In particular, vertical CP residuals between the case with and without GCPs highlight minimum values that similar to each other, while the maximum residual is slightly higher for the case without GCPs. Despite these considerations, the RMSE remains substantially the same (≈ 2 cm) for both configurations. All the values are reported in Table 3. Table 3. CP residuals for the oblique dataset (GCP residuals in brackets). In this case, the validation of DTMs has been made using 94 VPs, i.e., only those surveyed in the region covered by the double-grid flight plan. The standard deviation of the distribution of the elevation differences, even if slightly higher than those obtained for the nadiral dataset (≈ 5.5 vs. ≈ 4 cm), is comparable for the cases with and without the use of GCPs (5.3 vs. 5.5 cm, see Figure 9). The mean value also shows a difference between ∆H distributions that is less than 1 cm, so there seems to be no significant vertical offset for the oblique dataset between DTMs that consider GCPs and those that do not.

Nadiral Dataset with Oblique Images
The results of the alignment of the nadiral images in which a short sequence of oblique images has been added are reported in Table 4. The largest residuals can be found along the E direction, where they turn out to be slightly higher than those shown in Table 2 for the single-GCP configurations. Despite this, vertical residuals seem to be substantially unbiased, ranging from −0.061 to +0.043 m with an RMSE of 0.025 m. This value is half of the values obtained for the nadiral dataset without any oblique image in both single-GCP configurations (≈ 0.051 m).
The validation of the DTM, made using the same 112 VPs previously used for the validation of the nadiral dataset, highlights a distribution of vertical differences that is slightly worse than the case in which the nadiral dataset was processed introducing a set of GCPs. The mean value of the elevation differences ∆H changes from −2.3 to −3.2 cm, while the standard deviation of the distribution increases from 3.8 to 4.7 cm. However, the histograms reported in Figure 10 are certainly closer to the results obtained for the configuration with GCPs than the results of the single-GCP configurations (see Figure 8).

Cross Section Comparison
In geomorphological mapping, it is important to assess differences in cross-section profiles along the beach. For this reason, we extracted three profiles and then evaluated the bias between DTMs in two zones (see Figure 11). The first is the top of the dune, and the second is a part of the slope between the dune and the shoreline.
It is worth noting that the cross-sections have been selected in the region where nadiral and oblique flight plans overlapped, in order to compare all of the six DTMs generated in this work. Since we did not collect any data of cross-sections during survey operations, it is not possible to compare the profiles with a ground truth. However, the use of a set of GCPs produces a model with centimeter-level accuracy even without RTK/PPK camera locations (i.e., the traditional approach) that we assumed as reference. As regards the top of the dune, the most biased DTMs are those obtained with the introduction of a single GCP. The difference between them highlights an offset ranging from 0.08 to 0.10 m. Conversely, all the other DTMs have a lower difference, ranging from 0.01 to 0.05 m.
The second zone where we assessed the difference in DTMs confirms what happened at the top of the dune. The bias between nadiral single-GCP DTMs is up to 0.11 m, while other DTMs have a maximum difference of about 0.04 m.
The solution that integrated some oblique images into the nadiral dataset produced profiles that are generally unbiased and consistent with those generated using a set of GCPs, if we consider a vertical level of uncertainty of a few centimeters.

Discussion
The possibility of using UAVs without any GCPs is certainly a very interesting topic in many fields, as it would imply quicker image acquisition and the use of a limited number of instruments, i.e., only an aircraft using PPK or NRTK. The post-processing of RINEX files can bypass the need for a mobile data connection to receive NRTK corrections. Nonetheless, it makes the use of a base receiver benchmark unnecessary for RTK. However, the precision and accuracy of photogrammetric models depend on the quality of post-processed flight path solutions and decrease whenever the distance between the drone and the PPK station increases.
In this paper, we investigated the use of a CORS that is about 15 km away from our site. Such a distance is quite high for kinematic GNSS.
The results we obtained for the nadiral dataset show that RMSEs are particularly high along the E and U directions for single-GCP configurations, with CP RMSE values of 0.052 to 0.053 m and 0.051 to 0.052 m, respectively. In both cases, the 3D RMSE is about 0.075 m. One difference we noticed between the single-GCP models and that is further confirmed by the analysis made with VPs is the presence of a vertical offset between the models (and DTMs). The most probable reason is an error in estimating the focal length through the self-calibration in Agisoft Metashape, whose values for each instance of data processing are reported in Table 5. The configuration with a set of GCPs, indeed, is characterized by a focal length estimate that is approximately in the middle between the estimates of single-GCP cases. In addition, it shows a range of vertical residuals from a minimum of −0.037 m to a maximum of +0.038 m with an U RMSE of 0.021 m. The high quality of the DTM is also confirmed by VP analysis.
The results we obtained for the oblique dataset, in terms of CP residuals, highlight how the use of oblique images with PPK-based coordinates and the execution of double-grid missions provide excellent accuracy, even without any GCP. Both the maximum 3D residual (0.047 m) and the 3D RMSE (0.028 m) in such GCP-less configurations confirm the high quality of the photogrammetric model, also considering that the GSD is ≈ 2 cm. The focal length estimate for the oblique dataset is essentially the same whether a set of GCPs (3647.01 pix) is used or not (3646.74 pix). The validation of DTMs using the dataset of VPs shows mean values that are quite similar, even if standard deviations are double those of the RMSE obtained assessing the accuracy by CPs. However, it is worth noting that the actual duration of the VP survey through NRTK was higher than that of the signalization and survey of GCPs and CPs (see Section 2.1), which implies a broader dispersion of surveyed GNSS elevations for VPs. Finally, the introduction of a short sequence of oblique images within the nadiral dataset generated an alignment with good vertical accuracy (a 3D RMSE of 0.025 m), but the residuals in the E direction worsened (an E RMSE from ≈ 0.053 to 0.064 m). This produced an overall 3D RMSE that was quite high (0.070 m). The estimate of the focal length (3648.71 pix) is comparable to that obtained for the same dataset, without the use of a few oblique images and using a set of GCPs instead (3648.29 pix). The validation of the DTM of such a configuration highlighted a vertical difference that was greater than that of the oblique (and double-grid) dataset on average, in spite of a slightly lower standard deviation. However, the amount of VPs used for the validation is also lower (112 vs. 94).
The extraction of cross-section profiles from the DTMs and their further comparison has shown that most of the bias is present between models generated from nadiral images using a single GCP. The solution obtained with the introduction of a sequence of oblique images into the nadiral dataset, instead of a GCP, produced profiles that are generally comparable to those obtained using a set of GCPs. This is particularly clear for Cross Section B (see Figure 11), while more noise is present in Cross Section A.
Particularly interesting is the comparison of the results we reported here with those obtained by Taddia et al. [22] for the same datasets with RTK. In fact, RTK allows one to collect the camera center coordinates directly, at the time of each image capture and with a base-rover distance that is generally lower.
For what concerns the nadiral dataset, the cases in which a single GCP was used, either central or lateral, show an RMSE that for PPK is higher (3D RMSE = 0.075 m) than RTK (3D RMSE = 0.037 to 0.040 m). For the same dataset, the introduction of a set of well-distributed GCPs is able to compensate for the difference found between single-GCP configurations. 3D RMSEs evaluated for CPs become 0.029 m for RTK and 0.034 m for PPK. This result is a natural consequence of the traditional technique of using a set of GCPs to process images acquired without any RTK/PPK system. The subsequent validation of DTMs shows similar standard deviations of the distribution of vertical VP differences between the processed RTK and PPK introducing a single GCP (≈ 4 cm). However, it highlights a larger difference in terms of the mean value between cases of lateral and central GCPs. The poorer accuracy of PPK camera locations may affect this aspect, especially when the GCP is at the edge of the area and the photogrammetric block is "less stiff" because of the lower camera constraints. For both RTK and PPK, indeed, the mean values are −3.8 and −7 cm, respectively. Conversely, the use of a single GCP in a central position produced mean values of 0.1 (RTK) and 1.9 cm (PPK). Similar results between RTK and PPK were obtained in terms of VP elevation differences using a set of GCPs for the georeferencing of the models.
With regard to the oblique dataset, the comparison between RTK and PPK results shows very similar 3D RMSEs. GCP-less processing strategies are characterized by values of 0.030 and 0.028 m, respectively, while those using a set of GCPs are 0.027 and 0.025 m. The analysis along the E, N, and U directions highlights maximum differences between RTK and PPK RMSEs of about 1 to 1.5 cm.
Similarly, the comparison between DTMs generated using the oblique dataset with RTK and PPK shows a standard deviation (evaluated for the distribution of VP elevation differences) of about 5.5 cm for both RTK and PPK, with a mean value between −1.2 and +2 cm. This further confirms that DTMs by RTK and PPK oblique camera locations, with or without the use of GCPs, as well as the photogrammetric models from which they have been generated, are characterized by a similar level of accuracy.
The results show that DJI-P4RTK surveying methodologies exclusively based on the acquisition of aerial images (following a specific protocol, such as the recording of GNSS data, an accurate image capture timestamp, and a few oblique images in addition to nadiral ones) are able to produce photogrammetric models with a centimeter-level accuracy, even if no GCPs are placed or surveyed or if no local base receiver is used. This aspect becomes particularly relevant whenever the time that is actually available for the survey is very limited or the site is not accessible due to safety reasons [31,32], coastal pollution [33], or remote forested areas [34].

Conclusions
The results described in this work show how the DJI Phantom 4 RTK represents a tool that can, through the post-processing of the flight path GNSS data, provide all the elements necessary for the creation of photogrammetric models, DTMs, and cross-section profiles with centimeter-level accuracy.
The simple use of nadiral images acquired through single-grid missions with the introduction of a single GCP generates a DTM whose accuracy may be affected by a significant vertical offset. This aspect has been clearly highlighted by the results we obtained for the case in which we used a GCP at the edge of our site. On the other hand, the use of a set of GCPs eliminates the utility of having a multi-constellation, multi-frequency receiver on board.
The use of oblique images and the performance of double-grid flight plans ensure a higher level of precision without the need to use any GCPs. However, double-grid missions are certainly more time-consuming than single-grid ones. For this reason, we also investigated a "hybrid" situation in which we included a few oblique images within the nadiral dataset. The results show that the vertical accuracy of the photogrammetric model increased, even if it was not equal to the case of a complete oblique dataset acquired through double-grid missions.
It is also very interesting to note that there are no significant differences between the photogrammetric models and DTMs obtained for the oblique dataset without any GCP for RTK and PPK, even though the RAVE CORS is located about 15 km away from our site.
In the field of coastal mapping, the DJI-P4RTK certainly represents a surveying tool that is very promising and is able to provide centimeter-level accurate DTMs, even if a PPK strategy is followed. It is, therefore, easy to imagine that the fields of application of such a cheap and versatile drone will increase rapidly in the coming years or even months.  Acknowledgments: The authors are grateful to the student Alexandre Lazarou for his participation during GNSS surveying.

Conflicts of Interest:
The authors declare that there is no conflict of interest.