Spatial Structure, Short-temporal Variability, and Dynamical Features of Small River Plumes as Observed by Aerial Drones: Case Study of the Kodor and Bzyp River Plumes

Quadcopters can continuously observe ocean surface with high spatial resolution from relatively low altitude, albeit with certain limitations of their usage. Remote sensing from quadcopters provides unprecedented ability to study small river plumes formed in the coastal sea. The main goal of the current work is to describe structure and temporal variability of small river plumes on small spatial and temporal scales, which are limitedly covered by previous studies. We analyze optical imagery and video records acquired by quadcopters and accompanied by synchronous in situ measurements and satellite observations within the Kodor and Bzyp plumes, which are located in the northeastern part of the Black Sea. We describe extremely rapid response of these river plume to energetic rotating coastal eddies. We reveal several types of internal waves within these river plumes, measure their spatial and dynamical characteristics, and identify mechanisms of their generation. We suggest a new mechanism of formation of undulate fronts between small river plumes and ambient sea, which induces energetic lateral mixing across these fronts. The results reported in this study are addressed for the first time as previous related works were mainly limited by low spatial and/or temporal resolution of in situ measurements and satellite imagery.

General aspects of the structure and dynamics of river plumes as well as their regional features were addressed in many previous studies. Nevertheless, these works were mostly focused on large river plumes, while small rivers plumes received relatively little attention. However, small rivers play an important role in global land-ocean fluxes of fluvial water and suspended and dissolved sediments [54][55][56]. Small rivers form buoyant plumes that have small spatial scales and, therefore, small residence time of freshened water, which is equal to hours and days, due to relatively low volume of river discharge and its intense mixing with ambient sea [57]. Dissipation of freshened water as a result of mixing of a small plume with subjacent saline sea limitedly influences ambient sea and does not result in accumulation of freshwater in adjacent sea area. As a result, small plumes are characterized by sharp salinity and, therefore, density gradients at their boundaries with ambient sea. This feature is not typical for large river plumes and results in significant differences in spreading and mixing between small plumes and large plumes. Sharp vertical density gradient at the bottom boundary of a small plume hinders vertical energy transfer between a small river plume and subjacent sea [57]. This feature strongly affects spreading dynamics of a small plume due to the following reasons. First, the majority of wind energy transferred to sea remains in a small plume, because the vertical momentum flux diminishes at the density gradient between a plume and subjacent sea. Therefore, wind stress is concentrated in a shallow freshened surface layer that causes higher motion velocity and more quick response of dynamics of a small plume to variability of wind forcing, as compared to ambient sea [58,59]. It results in wind-driven dynamics of small plumes, which is characterized by very energetic short-temporal variability of their positions, shapes, and areas [60][61][62][63][64].
Study of structure and variability of small river plumes at small spatial and temporal scales is essential for understanding the fate of freshwater discharge from small rivers to sea and the related transport of suspended and dissolved river-borne constituents. However, high short-temporal variability of small plumes and their small vertical sizes inhibit precise in situ measurements of their thermohaline and dynamical characteristics [57]. Satellite remote sensing also does not provide the necessary spatial resolution and temporal coverage for small river plumes. As a result, many important aspects of structure, variability, and dynamics of small river plumes at small spatial and temporal scales remain unstudied.
Quadcopters are especially efficient in observation of small river plumes because they can continuously observe sea surface with high spatial resolution from relatively low altitude. Quadcopters can be used during overcast sky when optical satellite instruments cannot observe sea surface. The main drawback of their usage is relatively short duration of continuous operation (less than several hours), limited weight of carried instruments, and inability of their operation during strong wind, rain, snow, low temperature, and other inappropriate weather conditions. Despite these limitations, usage of quadcopters provides unprecedented ability to study structure of small river plumes, detect and measure their short-temporal variability, and register various dynamical features of these plumes. Therefore, the main goal of the current work is to describe structure and temporal variability of small river plumes on small spatial (from meters to hundreds of meters) and temporal (from minutes to hours) scales, which are limitedly covered by in situ measurements and satellite imagery and remain almost unaddressed by the previous studies.
In this work we use aerial remote sensing supported by synchronous in situ measurements and satellite observations to study small river plumes formed in the northeastern part of the Black Sea. We show that usage of aerial drones, first, strongly enhances in situ and satellite observations of structure and variability of small plumes, second, provides ability to perform accurate, continuous, and high-resolution measurements of their spatial characteristics and current velocity fields, and, finally, significantly improves operational organization of field measurements. Owing to continuous and high-resolution aerial remote sensing, we report several novel results about spatial structure, short-temporal variability, and dynamical features of small river plumes. These results include strongly inhomogeneous structures of small river plumes manifested by complex and dynamically active internal frontal zones; undulate (lobe-cleft) form of a sharp front between a small river plume and ambient sea; energetic lateral mixing across this front caused by its baroclinic instability; internal waves generated by river discharge near a river estuary and propagating within the inner plume; and internal waves generated by vortex circulation of river plume and propagating within the outer plume. The obtained results reveal significant differences in structure, variability, and dynamics between small plumes and large plumes.
The paper is organized as follows. Section 2 provides the detailed information about the aerial, in situ, and satellite data, as well as the processing methods used in this study. The results derived from aerial observations of small river plumes supported by in situ measurements and satellite observations are described in Section 3. Section 4 focuses on discussion and interpretation of the revealed features of spatial structure, short-temporal variability, and dynamics of small river plumes. The summary and conclusions are given in Section 5.

Study Area
In this work, we focused on the Kodor and Bzyp river plumes formed in the northeastern part of the Black Sea ( Figure 1). These rivers were chosen as the case sites due to the following reasons. First, these rivers have high concentrations of suspended sediments (300-500 g/m 3 in the Kodor River and 100-300 g/m 3 in the Bzyp River) [65], therefore, the turbid Kodor and Bzyp plumes can be effectively detected by optical aerial and satellite imagery. Second, the Kodor and Bzyp rivers are relatively small, their catchment areas are 2000 and 1500 km 2 , respectively, and their average annual discharges are approximately 130 and 120 m 3 /s, respectively [65]. As a result, the Kodor and Bzyp plumes are small enough to be observed by aerial remote sensing from relatively small altitude (< 200 m). However, both rivers are mountainous with large mean basin altitudes (> 1500 m) and slopes (> 0.02% ), as well as high drainage density (> 0.8 1/km). Therefore, during spring freshet and short-term rain-induced floods the runoffs from the Kodor and Bzyp rivers dramatically increase by 1-2 orders of magnitude. Third, despite their relatively small spatial extents, the Kodor and Bzyp plumes are the largest plumes in the study area. As a result, structure and dynamics of these plumes are not influenced by interaction with other river plumes. Fourth, the Kodor and Bzyp rivers have different mouth morphologies that affect the structure of their plumes. The majority of the Bzyp River runoff inflows to sea from the main river channel, however, a small side-channel is formed during high discharge periods. The Kodor River inflows to sea from three large river channels, which form the Kodor Delta. The mouths of these deltaic branches are located along the 2 km long segment of the coastline. Finally, wind, cloud, and rain conditions in the study area are favorable for aerial and satellite observations of the river plumes during the majority of the year.

Aerial, In Situ, and Satellite Data
Aerial observations of the Kodor and Bzyp plumes were performed by a quadcopter (DJI Phantom 4 Pro) equipped with a 12 MP/4K video camera. Aerial observations of the plumes were supported by ship-borne in situ measurements of salinity, temperature, turbidity, and current velocity within the plumes and the adjacent sea. The size of this quadcopter is small enough to be launched from and landed on a small boat. It provides opportunity for a quadcopter operator to be onboard the research vessel and to effectively coordinate synchronous in situ measurements and water sampling. Aerial observations and in situ measurements of the Kodor plume were conducted on 1-2 September 2018 and 1-3 April 2019, while aerial observations and in situ measurements of the Bzyp plume were performed on 31 May-1 June 2019. Below we provide the protocols of these aerial surveys according to the scheme suggested by Doukari et al. [74].
The quadcopter was flying over coastal sea areas adjacent to the Kodor and Bzyp river mouths.

Aerial, In Situ, and Satellite Data
Aerial observations of the Kodor and Bzyp plumes were performed by a quadcopter (DJI Phantom 4 Pro) equipped with a 12 MP/4K video camera. Aerial observations of the plumes were supported by ship-borne in situ measurements of salinity, temperature, turbidity, and current velocity within the plumes and the adjacent sea. The size of this quadcopter is small enough to be launched from and landed on a small boat. It provides opportunity for a quadcopter operator to be onboard the research vessel and to effectively coordinate synchronous in situ measurements and water sampling. Aerial observations and in situ measurements of the Kodor plume were conducted on 1-2 September 2018 and 1-3 April 2019, while aerial observations and in situ measurements of the Bzyp plume were performed on 31 May-1 June 2019. Below we provide the protocols of these aerial surveys according to the scheme suggested by Doukari et al. [74].
The quadcopter was flying over coastal sea areas adjacent to the Kodor and Bzyp river mouths. The take-off and landing spot was located on a vessel/boat that provided opportunity to perform flights at different areas of the plumes without any limitations on their distance to the seashore. The distance between the quadcopter and the research vessel/boat did not exceed 1 km. Quadcopter shooting altitude depended on the spatial scale of the sensing sea surface process and varied from 10-30 m for the small-scale frontal circulation to 150-200 m for detection of plume position and area. Weather conditions during the field surveys were favorable for usage of the quadcopter. Wind forcing during the flights was moderate (< 8 m/s), air temperature varied between 15 and 30 • C, and air humidity varied between 60% and 90%. The flights were conducted during no-rain conditions from morning to evening. In case of clear sky conditions, sun glint strongly affected the quality of the aerial data during the daytime. Wave heights were < 0.5 m during the flights.
In situ measurements performed in the study areas were the following. Continuous salinity and temperature measurements in the surface sea layer (0.5-1 m depth) were performed along the ship tracks using a shipboard pump-through system equipped by a conductivity-temperature-depth (CTD) instrument (Yellow Springs Instrument 6600 V2) [62,75]. Vertical measurements of salinity, temperature, and turbidity were performed using a CTD-instrument (Sea-Bird Electronics SBE 911plus) at 0.2 m spatial resolution. Vertical measurements of current velocity were performed using an acoustic Doppler current profiler (ADCP) (Teledyne RDI Workhorse Sentinel) and a CTD-ADCP-instrument (Aanderaa SeaGuard RCM). Vertical profiling was performed from sea surface to the depth of 10 m or to seafloor in shallow areas. The positions of individual in situ measurements are given in Section 3. Wind forcing during the field measurements was measured by a compact weather station (Gill GMX200) with temporal resolution of 1 minute. The weather station was mounted at the height of 10 m at a pier on a distance of 30 m from the coastline (Figure 1).

Processing of Aerial Data
In this study we used an optical flow algorithm to reconstruct velocity fields in the sea surface layer from quadcopter video records [76,77]. The main principle of optical flow algorithms used for calculation of motion from two consecutive pictures is the following. It was assumed that for each point → x (i.e., pixel) on both frames a certain signal intensity property I (i.e., brightness) was conserved: By linearizing the intensity of the second frame with respect to the intensity of the first frame a gradient constraint equation is obtained in the following way: where ∇I = I x , I y is the spatial partial derivatives of intensity, x /dt is the velocity, and I t is the temporal derivative of intensity. The derivatives ∇I and I t can be directly calculated, while the 2D velocity field → u is unknown. Therefore, Equation (1) requires an additional constraint and it is assumed Remote Sens. 2020, 12, 3079 6 of 30 that the displacement ∆ → x is constant in any small neighborhood, i.e., we search for a displacement that minimizes the constraint error: x is a weight function. Thus, minimization of E → x with respect to → u provides an additional condition for Equation (2). The resulting vector field → u calculated from Equations (2) and (3) is regarded as an optical flow estimate. In this work, we used the Farneback weight function [78] freely available in the OpenCV computer vision library (https://opencv.org/). This algorithm approximates a neighborhood of a pixel in each pair of frames by a quadratic polynomial function applying the polynomial expansion transform. Therefore, a constraint equation is based on a polynomial approximation of the given signal. On the assumption of small variability of a displacement field, the algorithm minimizes quadratic error of the constraint and calculates the optical flow estimation.
The estimation of surface velocity fields in the study region was performed in two stages. First, we applied the optical flow algorithm with large prescribed sizes of pixel neighborhoods for the reconstruction of motion of distinct plume boundaries and fronts. Second, we reconstructed motion within the river plume using the optical flow algorithm with a reduced neighborhood size. Spatial scale of motion, which is intended to be reconstructed, positively correlates with optimal size of a pixel neighborhood. An algorithm with a small pixel neighborhood more accurately reconstructs small-scale motion, but shows lower quality for large motion patterns, as compared to an algorithm with a large pixel neighborhood. The overall neighborhood size was prescribed according to spatial scales of ocean surface features (e.g., river plume fronts), in which motion is expected to be detected by an optical flow algorithm. Thus, the optimal neighborhood size intended to reconstruct the large-scale motion of river plumes should be equal to the width of distinct plume boundaries and fronts. In this study, the large size of a pixel neighborhood was prescribed equal to 30 m, while the small size of a pixel neighborhood was set equal to 1 m. In case of application of this algorithm to other regions, we suggest prescribing neighborhood sizes equal to relevant spatial scales of the considered ocean surface features.
Due to high resolution of the video camera used and continuous video recording, the optical flow algorithm efficiently detected motion of the distinct frontal zones within the river plumes, as well as motion of foam and floating litter accumulated at these fronts which is indicative of the circulation patterns at the frontal zones. As a result, the reconstructed surface velocity fields showed good accordance with visually inspected shifts of the frontal zones, foam, and floating litter at the video records. Stable positioning of a quadcopter is important for precise motion detection at sea surface. Moderate wind speed during the field surveys did not negatively affect the quality of the obtained aerial data. However, strong wind forcing during camera shooting can hinder accurate reconstruction of surface velocity fields. Sun glint is another important issue that can impede motion detection at aerial video records. Intensity of the sun glint depends on solar elevation angle, camera shooting angle and direction; therefore, it can be reduced by correct selection of quadcopter altitude and position. Usage of polarizing filters for quadcopter camera can reduce glint from water surface, however, its efficiency strongly depends on camera shooting angle.

Spatial Structure and Short-temporal Variability of the Kodor and Bzyp Plumes
The field surveys were performed during spring freshet at the Bzyp River Vertical salinity measurements in the study areas revealed that these low-saline plumes are shallow (< 5 m depth) and have distinct vertical salinity gradients with the ambient saline sea. Due to elevated concentrations of terrigenous suspended sediments in the Kodor and Bzyp rivers [65], turbidity within the Kodor and Bzyp plumes was significantly larger than in the ambient sea and showed good correlation with reduced salinity (Figure 2). The Pearson correlation coefficients (r) between salinity and turbidity are equal to −0.87 and −0.71 for the Kodor and Bzyp plumes respectively with p-values equal to 0.0000. These high absolute values of the correlation coefficients at low p-values indicate that the observed relations between salinity and turbidity within the Kodor and Bzyp plumes (low salinity and high turbidity), on the one hand, and the ambient sea water (high salinity and low turbidity), on the other hand, are statistically significant. As a result, surface turbidity structures of the Kodor and Bzyp plumes observed by optical remote sensing are indicative of surface salinity structures of these plumes.   Aerial remote sensing and satellite imagery showed that the alongshore extents of turbid surface water associated with the considered river plumes during low discharge conditions are 1−5 km. The obtained estimates were consistent with salinity measurements at the study area. However, flooding discharge results in abrupt expanding of these plumes, their extents and areas can exceed 20 km and 50 km 2 , respectively. Aerial and satellite images, surface salinity distribution, and vertical salinity profiles obtained on 31 August 2018 in the coastal area adjacent to the Kodor Delta are illustrative of spatial scales, as well as horizontal and vertical structure of the Kodor plume ( Figure 3).
Aerial observations and in situ measurements revealed strongly inhomogeneous salinity and turbidity structure of the Kodor plume manifested by complex and dynamically active frontal zones within the plume (Figures 4-6). In particular, surface salinity showed no dependence on the distance to the mouths of the deltaic branches that is regarded typical for river plumes [79][80][81], especially in numerical modeling studies [82][83][84][85]. This inhomogeneous structure is formed due to impact of several different processes including the formation of the Kodor plume by several spatially distributed sources, the large inter-day river discharge variability in response to sporadic rain events, and the bathymetric features that influence spreading of the plume.  The Kodor River inflows to sea from three deltaic branches with different discharge rates. As a result, all three branches form individual river plumes that merge and coalesce into the common Kodor plume. These three river plumes have different structure, spatial characteristics, and dynamics, therefore, they interact as individual water masses and form stable frontal zones observed by aerial imagery (Figure 4a) [86][87][88]. In situ measurements performed on 2 September, 2018 revealed sharp salinity gradient at the frontal zone between the river plumes formed by the northern and the central deltaic branches of the Kodor River. Surface salinity along the transect that crossed this frontal zone abruptly decreased from 14 to 8-10 on a distance of 5 m (Figure 4b).
turbidity structure of the Kodor plume manifested by complex and dynamically active frontal zones within the plume (Figures 4-6). In particular, surface salinity showed no dependence on the distance to the mouths of the deltaic branches that is regarded typical for river plumes [79][80][81], especially in numerical modeling studies [82][83][84][85]. This inhomogeneous structure is formed due to impact of several different processes including the formation of the Kodor plume by several spatially distributed sources, the large inter-day river discharge variability in response to sporadic rain events, and the bathymetric features that influence spreading of the plume. The Kodor River inflows to sea from three deltaic branches with different discharge rates. As a result, all three branches form individual river plumes that merge and coalesce into the common Kodor plume. These three river plumes have different structure, spatial characteristics, and The discharge of the Kodor River shows quick response to precipitation events that is common for small mountainous rivers with small and steep watershed basins. Frequent rains at the mountainous northeastern coast of the Black Sea cause high inter-day variability of the discharge rate of the Kodor River [65,89]. As a result, the area of the Kodor plume can significantly change during less than one day that was observed on 31 August-2 September 2018 during the field survey. Heavy rain that occurred during 6 hours at night on 31 August-1 September (according to the local weather station measurements) caused increase of the river discharge from 80 to 150 m 3 /s during several hours. The area of the Kodor plume doubled from 31 August to 1 September in response to the flash flood. Wind direction during 31 August-1 September was stable (southwestern), while wind speed slightly increased from 2-3 m/s to 4-5 m/s. Then river discharge steadily decreased to pre-flooding conditions, which were registered on 2 September, while wind direction changed to eastern and wind velocity decreased to 3-4 m/s. In situ measurements and aerial remote sensing performed on 2 September, i.e., shortly after the flood, observed, first, the large residual plume that was formed on 1 September during the flooding event and did not dissipate yet and, second, the emergent plume that was formed on 2 September after the decrease of river discharge rate ( Figure 5). These plumes had different spatial scales, structures, thermohaline, and dynamical characteristics. As a result, similarly to the river plumes formed by different deltaic branches, the residual and the emergent plumes interacted as individual water masses and formed complex frontal zones within the common Kodor plume. performed on 2 September, i.e., shortly after the flood, observed, first, the large residual plume that was formed on 1 September during the flooding event and did not dissipate yet and, second, the emergent plume that was formed on 2 September after the decrease of river discharge rate ( Figure 5). These plumes had different spatial scales, structures, thermohaline, and dynamical characteristics. As a result, similarly to the river plumes formed by different deltaic branches, the residual and the emergent plumes interacted as individual water masses and formed complex frontal zones within the common Kodor plume. Interaction between the Kodor plume and the seafloor at the shallow zones is the third process that induces inhomogeneous structure of this plume. Aerial imagery detected the area of reduced turbidity formed behind the shoal, which is located in front of the northern deltaic branch ( Figure 6). This low-turbid zone contrasted especially with the surrounding turbid river plume during the flooding discharge on 1 September 2018. In situ measurements showed that surface salinity at this low-turbid zone (15) was significantly greater than at the adjacent turbid part of the plume (12.5-13) (Figure 6c). Surface circulation also differed in these two parts of the plume. The northward flow (10 cm/s) was observed in the low-turbid zone, while the southeastward flow (20 cm/s) dominated in the adjacent turbid part (Figure 6d). The formation of this zone is caused by the interaction of the inflowing river jet with seafloor at the shoal that induces deceleration of the jet and its increased mixing with saline and low-turbid sea water. The stable front bounding this low-turbid and high-saline zone inside the plume was observed on a distance of up to 1 km from the shoal.
low-turbid zone (15) was significantly greater than at the adjacent turbid part of the plume (12.5-13) (Figure 6c). Surface circulation also differed in these two parts of the plume. The northward flow (10 cm/s) was observed in the low-turbid zone, while the southeastward flow (20 cm/s) dominated in the adjacent turbid part (Figure 6d). The formation of this zone is caused by the interaction of the inflowing river jet with seafloor at the shoal that induces deceleration of the jet and its increased mixing with saline and low-turbid sea water. The stable front bounding this low-turbid and high-saline zone inside the plume was observed on a distance of up to 1 km from the shoal.

Dynamical Features of the Kodor and Bzyp Plumes
Using aerial remote sensing we detected several dynamical features of the Kodor and Bzyp plumes and measured their spatial characteristics. Based on the surface velocity data reconstructed from the aerial video records, we studied dynamical characteristics of these features and analyzed

Dynamical Features of the Kodor and Bzyp Plumes
Using aerial remote sensing we detected several dynamical features of the Kodor and Bzyp plumes and measured their spatial characteristics. Based on the surface velocity data reconstructed from the aerial video records, we studied dynamical characteristics of these features and analyzed their physical background. Aerial remote sensing detected a swirling eddy within the Kodor plume on 1 September 2018 (Figure 7). This eddy was formed at the southern part of the emergent plume at its border with the residual plume near the Iskuria Cape. The aerial image of this part of the plume acquired at 12:52 ( Figure 7a) showed inhomogeneous structure of the emergent plume without any eddy. The distinct border between the emergent and the residual plumes was stretched from the Iskuria Cape in the northwestern direction. The beginning of formation of the eddy was registered at 14:42 (Figure 7b), then at 15:34 the well-developed eddy was observed (Figure 7c,d). The diameter of the eddy was approximately 500 m, it was rotating in an anticyclonic direction, while its center was moving at an angle of approximately 30 • across the border of the emergent plume. Processing of the video record of this eddy provided estimations of velocity of its movement (0.9 m/s) and rotation (0.4 m/s). The aerial observations performed at 16:16 did not show any surface manifestations of the eddy at the study area; therefore, we presume that it shifted off the observation area during less than an hour. Wind conditions were stable during the considered period, wind speed did not exceed 3.5 m/s.  In situ thermohaline and velocity measurements were performed within the eddy at 15:57-16:01 (Figure 7e,f). They included continuous measurements at a depth of 0.7-0.8 m for 4.5 minutes followed by vertical profiling from surface to the depth of 13 m. Note that the measurements were performed at the stable point, while the eddy was moving. As a result, the performed measurements registered salinity and velocity in different parts of the eddy while it was passing the point of measurements. The intense northward flow (55 cm/s) registered in the surface layer at the beginning of the measurements steadily dissipated to <10 cm/s during the first stage of the measurements (Figure 7f). The eastward velocity component was slightly positive during the first two minutes of the measurements (6 cm/s on average with the peak value of 16 cm/s) and then changed to slightly negative (−5 cm/s on average with the peak value of −11 cm/s). It was accompanied by significant variability of salinity that increased from 13.5 to 15.5 during the first 1.5 min of the measurements and then decreased to 13.5 (Figure 7e). The observed variability of velocity and salinity in the surface layer confirms northward propagation and anticyclonic rotation of this eddy observed at aerial video (Supplementary Materials). However, the movement and rotation velocities registered by in situ measurements were twice less than those reconstructed from the aerial video. This difference is caused by the fact that in situ measurements were performed not at the central part of the eddy, but at its periphery. The observed variability of salinity in the surface layer was caused by intrusion of saline water from the ambient sea to the plume induced by the rotation of the eddy (Figure 7d). Vertical profiles of salinity and velocity measured at 16:02, i.e., after the measurements in the surface layer, registered strong northwestward flow in the subjacent saline sea (Figure 7g,h). Its maximal velocity (15-25 cm/s) was observed immediately beneath the plume at depths of 3-5 m, then velocity decreased to 10-15 cm/s at depths of 8-9 m and to <5 cm/s at depths of 10-13 m. This northwestward flow (20-30 cm/s) was also registered along the Iskuria Cape at the previous day that confirms the presence of the northwestward jet behind the Iskuria Cape which is presumed to generate the observed eddy.
Interaction between sub-mesoscale eddies and the Kodor plume was also observed by satellite imagery. The chains of small anticyclonic eddies (300-500 m in diameter) formed behind the Iskuria Cape and interacting with the Kodor plume were registered on 17 July 2018, 21 August 2019, and 26 August 2019 (Figure 8a). Positions, sizes, and shapes of four to five subsequent eddies within these chains indicate that these chains were periodically generated near the Iskuria Cape and propagated in the northwestward direction shortly before the periods of satellite observations. While tracks of the eddies were crossing the Kodor plume, the turbid plume water was twisted into the eddies, which made them visible at satellite imagery. After these eddies propagated off the plume the trapped turbid water remained connected with the plume that illustrated difference in trajectories and velocities of the eddies and the wind-driven far-field part of the plume (Figure 8a).
Satellite images acquired during the periods of field measurements at the Kodor plume did not register interactions between the eddies and the plume due to episodic character of these features, i.e., eddies do not constantly form and propagate at the study area. Therefore, the satellite images presented in Figure 8 are not synchronous with the field surveys. However, sizes and anticyclonic rotation in the northwestward direction were similar for eddies detected at the Kodor plume by aerial and satellite remote sensing. As a result, we presume that we observe the same process and, therefore, can jointly analyze its spatial and temporal characteristics obtained from aerial and satellite measurements. Satellite imagery also observed eddies formed behind the Pitsunda Cape and interacting with the Bzyp plume on 30 July 2017 and 10 October 2019 (Figure 8b). However, in contrast to the eddies registered within the Kodor plume, these eddies were individual, i.e., did not form chains. Moreover, these eddies were much larger (2-4 km in diameter) and were rotating in cyclonic direction. Satellite images acquired during the periods of field measurement at the Bzyp plume also did not register interactions between eddies with the Bzyp plume. Satellite images acquired during the periods of field measurements at the Kodor plume did not register interactions between the eddies and the plume due to episodic character of these features, i.e., eddies do not constantly form and propagate at the study area. Therefore, the satellite images presented in Figure 8 are not synchronous with the field surveys. However, sizes and anticyclonic rotation in the northwestward direction were similar for eddies detected at the Kodor plume by aerial and satellite remote sensing. As a result, we presume that we observe the same process and, therefore, can jointly analyze its spatial and temporal characteristics obtained from aerial and satellite measurements. Satellite imagery also observed eddies formed behind the Pitsunda Cape and interacting with the Bzyp plume on 30 July 2017 and 10 October 2019 (Figure 8b). However, in contrast to the eddies registered within the Kodor plume, these eddies were individual, i.e., did not form chains. Moreover, these eddies were much larger (2-4 km in diameter) and were rotating in cyclonic direction. Satellite images acquired during the periods of field measurement at the Bzyp plume also did not register interactions between eddies with the Bzyp plume.
Satellite image acquired on 10 October 2019 detected packets of internal waves emerging from the rotating eddy and propagating within the Bzyp plume (Figure 9b). Aerial observations on 1  (Figure 9b) are not synchronized and show different river plumes at different dates. Aerial and satellite images acquired during the period of field measurements at the Bzyp plume did not register internal waves within the Bzyp plume. Therefore, in Figure 9 we show airborne images of internal waves at the Kodor plume and satellite images of internal waves at the Bzyp plume.
were generated within the plume around the rotating eddy. On the other hand, airborne remote sensing provided opportunity to detect individual internal waves with high spatial resolution and to register their velocities (Figure 9a). High-resolution aerial imagery detected that the distances between the individual waves within the wave packet in the Kodor plume were 2-4 m. The length of the wave packet front was approximately 200 m. The number of waves within the wave packet varied from 12 at its northern part to 3 at its southern periphery. Processing of high-resolution video records revealed that velocity of the wave packet was equal to 0.21 m/s.  Osadchiev [33] described a mechanism of generation of internal waves in small river plumes as a result of rapid deceleration of an inflowing river jet and formation of a hydraulic jump in vicinity of a river mouth. These internal waves propagate offshore and are regularly observed by satellite imagery in many coastal regions in the World [33,90,91]. Using aerial remote sensing we recorded generation and propagation of these internal waves from the mouth of the side-channel of the Bzyp River on 1 July 2019 (Figure 11a). The internal waves were generated at a distance of 40-50 m from the river mouth every 19 seconds on average, i.e., 29 individual waves were generated during a 9-min long video recording of this area. The distances between the waves decreased from 8-10 m near the river mouth to 1-2 m at the distance of 500 m from the river mouth. Wave velocities were equal to 0.27-0.31 m/s. Moderate (2-3 m/s) northern wind was registered during the considered period.
Aerial observations of internal waves in the Bzyp plume described above were supported by in situ salinity and turbidity measurements performed from a flat-bottomed boat with shallow draft to minimize the boat-induced mixing of sea surface layer (Figure 11). Measurements included 15 surface-to-bottom profiles continuously performed from a free-drifting boat starting at the Despite a large difference in coverage and spatial resolution of the aerial and satellite imagery presented in Figure 9, they both distinctly demonstrate propagation of internal waves within the river plumes. Satellite remote sensing has wide spatial coverage and provides information about spatial characteristics of wave packets at different parts of the plumes (Figure 9b). Distances between the wave packets observed at Sentinel-2 satellite images varied from 30 to 150-200 m, while lengths of the wave packets were up to 5-6 km. Satellite images demonstrated that dozens of internal waves were generated within the plume around the rotating eddy. On the other hand, airborne remote sensing provided opportunity to detect individual internal waves with high spatial resolution and to register their velocities (Figure 9a). High-resolution aerial imagery detected that the distances between the individual waves within the wave packet in the Kodor plume were 2-4 m. The length of the wave packet front was approximately 200 m. The number of waves within the wave packet varied from 12 at its northern part to 3 at its southern periphery. Processing of high-resolution video records revealed that velocity of the wave packet was equal to 0.21 m/s. Aerial remote sensing also detected multiple packets of low-frequency internal waves that propagated within the Kodor plume towards the coast on 2 September 2018 ( Figure 10). These packets consisted of 5-15 waves that were stretched along the coast, albeit had complex shapes not related to the shapes of the plume front or the coastline. Distances between individual waves varied from 5 to 70 m in the observed wave packets. Frontal length of these packets varied from~100 m (Figure 10a) to 2-3 km (Figure 10b), while their speeds were 10-15 cm/s. Wind speed during this period was 2-3 m/s. Osadchiev [33] described a mechanism of generation of internal waves in small river plumes as a result of rapid deceleration of an inflowing river jet and formation of a hydraulic jump in vicinity of a river mouth. These internal waves propagate offshore and are regularly observed by satellite imagery in many coastal regions in the World [33,90,91]. Using aerial remote sensing we recorded generation and propagation of these internal waves from the mouth of the side-channel of the Bzyp River on 1 July 2019 (Figure 11a). The internal waves were generated at a distance of 40-50 m from the river mouth every 19 seconds on average, i.e., 29 individual waves were generated during a 9-min long video recording of this area. The distances between the waves decreased from 8-10 m near the river mouth to 1-2 m at the distance of 500 m from the river mouth. Wave velocities were equal to 0.27-0.31 m/s. Moderate (2-3 m/s) northern wind was registered during the considered period.
Aerial observations of internal waves in the Bzyp plume described above were supported by in situ salinity and turbidity measurements performed from a flat-bottomed boat with shallow draft to minimize the boat-induced mixing of sea surface layer ( Figure 11). Measurements included 15 surface-to-bottom profiles continuously performed from a free-drifting boat starting at the generation area of the internal waves at the distance of 10 m from the river mouth and finishing 90 m far from the starting point (Figure 11a). The obtained data revealed large difference in vertical salinity structure of the Bzyp plume inside and outside this generation area of internal waves. The first half of the hydrological transect was located at the area of formation of the hydraulic jump as a result of abrupt deceleration of the inflowing river jet (Figure 11b). Similarly to the hydraulic jump observed and described by Osadchiev [33] at the inflowing jet of the Mzymta River, we registered anomalously deep penetration of low-saline water at the generation area of the internal waves in the Bzyp plume. Low-saline water (10)(11)(12)(13)(14) was observed from surface to the depth of 3-4 m along 0-5 m and 25-35 m of the transect. Vertical salinity structure within this part of the plume was unstable with multiple overturns (reverse salinity difference was up to 1 at vertical distance of 0.1 m) and large salinity gradients. Vertical salinity structure of the Bzyp plume between the areas of the hydraulic jumps, i.e., along the 5-25 m of the transect, showed relatively homogenous salinity (14.5-16) from surface to bottom, albeit it was much higher than within the areas of hydraulic jumps.
Outside the generation area of the internal waves, i.e., along the 35-90 m of the transect, surface salinity was relatively homogenous (14.5-15.5) and vertical salinity structure was stable. Vertical salinity gradient outside the generation area of internal waves was two orders of magnitude less than the largest values registered in the hydraulic jumps. However, salinity measurements did not cover top 0.5 m of the surface layer, where presumably was located the salinity gradient. Vertical turbidity structure, however, did not show large difference inside and outside the generation area of the internal waves (Figure 11c). The turbid layer was observed from surface to the depth of 1-1.5 m along the first part of the transect and then its depth steadily decreased to 0.5 m. This feature shows that salinity and turbidity structure of a river plume can be significantly different in areas of very intense advection and turbulent mixing. the internal waves (Figure 11c). The turbid layer was observed from surface to the depth of 1-1.5 m along the first part of the transect and then its depth steadily decreased to 0.5 m. This feature shows that salinity and turbidity structure of a river plume can be significantly different in areas of very intense advection and turbulent mixing.

Undulate Borders of the Kodor and Bzyp plumes
Aerial remote sensing of the Kodor and Bzyp plumes showed undulate structure of long segments of their outer borders manifested by alternation of specific convex and concave segments. These segments are 2-10 m long and up to 2 m wide and hereafter are referred as "lobes" and "clefts" [52,53]. Aerial images of the undulate fronts observed at the Kodor plume border on 1 September, 2018 and at the Bzyp plume border on 1 June 2019 are shown in Figure 12. This lobe-cleft structure was registered only at sharp and narrow frontal zones formed between the emerging plume, on the one hand, and the residual plume or the ambient sea, on the other hand. Lobes and clefts were absent at diffuse fronts, i.e., wide and low-gradient fronts that contour the outer parts of the plumes, which experience intense mixing with the ambient sea. In particular, these undulate fronts commonly extended from the river mouths and bounded the inflowing river jets, i.e., near-field parts of the plumes. These fronts were not observed in the far-field parts of the plumes and in the coastal surf zone during periods of active wave breaking due to intense mixing ( Figure 12).
Remote Sens. 2020, 12, x FOR PEER REVIEW 17 of 29 Figure 11. (a) aerial image of surface manifestations of internal waves propagating within the Bzyp plume off the river mouth and location of the hydrological transect (black line) on 1 July 2019 (a). Vertical (b) salinity and (c) turbidity profiles along the hydrological transect.

Undulate Borders of the Kodor and Bzyp plumes
Aerial remote sensing of the Kodor and Bzyp plumes showed undulate structure of long segments of their outer borders manifested by alternation of specific convex and concave segments. These segments are 2-10 m long and up to 2 m wide and hereafter are referred as "lobes" and "clefts" [52,53]. Aerial images of the undulate fronts observed at the Kodor plume border on 1 September, 2018 and at the Bzyp plume border on 1 June 2019 are shown in Figure 12. This lobe-cleft structure was registered only at sharp and narrow frontal zones formed between the emerging plume, on the one hand, and the residual plume or the ambient sea, on the other hand. Lobes and clefts were absent at diffuse fronts, i.e., wide and low-gradient fronts that contour the outer parts of the plumes, which experience intense mixing with the ambient sea. In particular, these undulate fronts commonly extended from the river mouths and bounded the inflowing river jets, i.e., near-field parts of the plumes. These fronts were not observed in the far-field parts of the plumes and in the coastal surf zone during periods of active wave breaking due to intense mixing ( Figure  12). We observed significant short-temporal variability of the undulate fronts induced by the following recurrent process (Figure 13). Once a lobe is formed, it starts to increase seaward. Ballooning of neighboring lobes results in their coalescence and the subsequent merging. At the same time the cleft between these lobes is steadily decreasing and transforms into a spot of saline sea (with area of 0.1-0.5 m 2 ) isolated from the ambient sea, i.e., trapped by the merged lobes within the plume (Figure 13). The merged lobes and the trapped saline sea area finally dissipate, and then the process of formation of new lobes at this part of the plume front restarts. The continuous recurrent process of formation of lobes, their merging, and subsequent dissipationwas observed along the undulate fronts of the Kodor and Bzyp plumes. Residual time of an individual lobe, i.e., from its formation to dissipation, was 1-2 min. We observed significant short-temporal variability of the undulate fronts induced by the following recurrent process (Figure 13). Once a lobe is formed, it starts to increase seaward. Ballooning of neighboring lobes results in their coalescence and the subsequent merging. At the same time the cleft between these lobes is steadily decreasing and transforms into a spot of saline sea (with area of 0.1-0.5 m 2 ) isolated from the ambient sea, i.e., trapped by the merged lobes within the plume ( Figure 13). The merged lobes and the trapped saline sea area finally dissipate, and then the process of formation of new lobes at this part of the plume front restarts. The continuous recurrent process of formation of lobes, their merging, and subsequent dissipationwas observed along the undulate fronts of the Kodor and Bzyp plumes. Residual time of an individual lobe, i.e., from its formation to dissipation, was 1-2 min.  Due to convergence of surface currents at sharp plume fronts [92], foam and floating litter commonly accumulate at the undulate fronts of the plumes (Figures 13a and 14a). Using optical flow processing of aerial video records, we detected motion of foam and floating litter and reconstructed surface circulation along the undulate fronts of the Kodor and Bzyp plumes ( Figure 14). The circulation structure within the lobes consists of pairs of cyclonic and anticyclonic vortices that form, balloon, merge, and dissipate with the lobes (black lines in Figure 14b). The trajectories of foam and floating litter revealed that cyclonic vortices are significantly more prominent and intense, as compared to anticyclonic vortices. Foam and floating litter are mainly accumulated within cyclonic eddies, i.e., in the right parts of the lobes if we look from the sea towards the plume (Figure 14a). Foam and floating litter are rotated by cyclonic eddies within the right parts of the lobes during the majority of time of aerial observations. Once a parcel of foam or floating litter is advected off a cyclonic eddy and enters an anticyclonic eddy in the left part of the lobe, it is transported to the outer part of the lobe and then is trapped by the cyclonic eddy in the neighboring (leftward) lobe (red lines in Figure 14b). As a result, these parcels are skipping leftward between the right parts of lobes. Therefore, foam and floating litter are steadily transported to the left along the plume border. The observed large intensity of cyclonic circulation within the lobes, as compared to anticyclonic circulation, is presumed to have the same background as the dominance of cyclonic spirals at satellite images of sea surface caused by differences between the rotary characteristics of cyclonic and anticyclonic eddies in the sea [93].
Remote Sens. 2020, 12, x FOR PEER REVIEW 19 of 29 Due to convergence of surface currents at sharp plume fronts [92], foam and floating litter commonly accumulate at the undulate fronts of the plumes (Figures 13a and 14a). Using optical flow processing of aerial video records, we detected motion of foam and floating litter and reconstructed surface circulation along the undulate fronts of the Kodor and Bzyp plumes ( Figure 14). The circulation structure within the lobes consists of pairs of cyclonic and anticyclonic vortices that form, balloon, merge, and dissipate with the lobes (black lines in Figure 14b). The trajectories of foam and floating litter revealed that cyclonic vortices are significantly more prominent and intense, as compared to anticyclonic vortices. Foam and floating litter are mainly accumulated within cyclonic eddies, i.e., in the right parts of the lobes if we look from the sea towards the plume (Figure 14a). Foam and floating litter are rotated by cyclonic eddies within the right parts of the lobes during the majority of time of aerial observations. Once a parcel of foam or floating litter is advected off a cyclonic eddy and enters an anticyclonic eddy in the left part of the lobe, it is transported to the outer part of the lobe and then is trapped by the cyclonic eddy in the neighboring (leftward) lobe (red lines in Figure 14b). As a result, these parcels are skipping leftward between the right parts of lobes. Therefore, foam and floating litter are steadily transported to the left along the plume border. The observed large intensity of cyclonic circulation within the lobes, as compared to anticyclonic circulation, is presumed to have the same background as the dominance of cyclonic spirals at satellite images of sea surface caused by differences between the rotary characteristics of cyclonic and anticyclonic eddies in the sea [93].  We presume that the undulate structure of the sharp plume borders is formed due to baroclinic instability between the plumes and the ambient sea. The pressure gradient force across the front is equal to where g is the gravity acceleration, ∆ρ is the density difference between the plume and the ambient sea, ρ sea is the density of the sea, h is the depth of the plume, and x is the cross-front direction. In situ measurements performed at the undulate fronts showed that surface salinity abruptly increased across these fronts ( Kodor plume at the narrow frontal zone was 2 m (Figure 15b), the depth of the Bzyp plume was 4 m. As a result, the values of pressure gradient across these frontal zones calculated from Equation (4) are equal to 0.05 and 0.1 m/s 2 for the Kodor and Bzyp plumes, respectively. This large pressure gradient observed across the plume fronts is the source of potential energy that induces formation of lobes and clefts as follows. Small perturbation of a sharp frontal zone and the subsequent formation of a local convex segment cause increase of local length of the front and, therefore, increase of the cross-front advection induced by the pressure gradient. It results in ballooning of the lobe till it coalesces and merges with the neighboring lobe. Merging of two lobes accompanied by trapping of a spot of saline sea water and its subsequent mixing with the plume water cause a reduction of local salinity anomaly and, therefore, a decrease of local pressure gradient. It hinders formation of a lobe at this segment of the plume, while new lobes are formed at the adjacent segments of the plume front. Therefore, baroclinic instability causes formation, merging, and dissipation of the observed lobe-cleft structures and influences mixing between the river plumes and the ambient sea.
Aerial imagery detected the 3-4 m wide stripe of low-turbid water within the Kodor plume located at the distance of 10-20 m from the undulate border and stretched along this border ( Figure  15a). We presume that this low-turbid stripe is formed as a result of continuous trapping of spots of saline sea water by merging lobes. Horner-Devine et al. [53] assumed that the lobe-cleft structure is formed by subsurface vortexes that are propagating from the inner part of the plume towards its border with the ambient sea. However, aerial video records showed stable position and shape of this stripe that evidences absence of any subsurface vortexes described by Horner-Devine et al. [53].  This large pressure gradient observed across the plume fronts is the source of potential energy that induces formation of lobes and clefts as follows. Small perturbation of a sharp frontal zone and the subsequent formation of a local convex segment cause increase of local length of the front and, therefore, increase of the cross-front advection induced by the pressure gradient. It results in ballooning of the lobe till it coalesces and merges with the neighboring lobe. Merging of two lobes accompanied by trapping of a spot of saline sea water and its subsequent mixing with the plume water cause a reduction of local salinity anomaly and, therefore, a decrease of local pressure gradient. It hinders formation of a lobe at this segment of the plume, while new lobes are formed at the adjacent segments of the plume front. Therefore, baroclinic instability causes formation, merging, and dissipation of the observed lobe-cleft structures and influences mixing between the river plumes and the ambient sea.
Aerial imagery detected the 3-4 m wide stripe of low-turbid water within the Kodor plume located at the distance of 10-20 m from the undulate border and stretched along this border (Figure 15a). We presume that this low-turbid stripe is formed as a result of continuous trapping of spots of saline sea water by merging lobes. Horner-Devine et al. [53] assumed that the lobe-cleft structure is formed by subsurface vortexes that are propagating from the inner part of the plume towards its border with the ambient sea. However, aerial video records showed stable position and shape of this stripe that evidences absence of any subsurface vortexes described by Horner-Devine et al. [53].

Discussion
In this study, we obtained several important results about structure, short-temporal variability, and dynamics of small river plumes. First, we revealed strongly inhomogeneous structure of small plumes manifested by multiple frontal zones between different parts of the plumes. These parts have different structures and dynamical characteristics and interact as individual water masses. Second, we reported fast motion of small plumes caused by interaction with coastal eddies. Third, we observed generation and propagation of different types of internal waves within small plumes. Forth, we described formation of lobe-cleft structures at sharp borders of small plumes and reported intense lateral mixing across these fronts caused by their baroclinic instability. The results listed above are important for understanding spreading and mixing of small plumes, however, they are addressed for the first time as previous related works were mainly limited by low spatial and/or temporal resolution of in situ measurements and satellite imagery. Below we provide physical interpretation of these features observed at the Kodor and Bzyp plumes and discuss importance of their study at other small plumes in the World Ocean.
In general, river plumes are regarded as "smooth" water masses without internal fronts and sharp gradients. This approach is widely used in analytical and numerical modeling studies focused on river plumes, including the fundamental and highly cited papers [82][83][84][85][94][95][96]. Many relevant studies based on in situ and satellite data confirmed that this approach provides realistic results for buoyant plumes formed by large rivers plumes which internal structures indeed are characterized by steady changes of salinity and other characteristics. In this work, we present the results of aerial remote sensing of the Kodor and Bzyp plumes supported by in situ measurements that provide an evidence of strongly inhomogeneous internal structure of small plumes. This structure is manifested by complex internal frontal zones and sharp salinity and turbidity gradients within small plumes. These gradients and frontal zones strongly modify circulation within the plumes, in particular, they hinder cross-frontal advection within the plumes and separate them into semi-isolated, but interacting structures. Therefore, identification and study of the processes that govern formation of frontal zones within small plumes is important for understanding of spreading and mixing of freshwater discharge in the sea and the related transport of river-borne suspended and dissolved material.
The Kodor River inflows to the Black Sea from multiple deltaic branches and forms several river plumes. These plumes are closely located; they interact as individual water masses and coalesce into the common Kodor plume. Interaction, collision, and coalescence of buoyant plumes formed by rivers, which estuaries are located in close proximity, were addressed in several previous studies [86][87][88][97][98][99]. Similar processes occur within plumes formed by freshwater discharge from multiple deltaic branches, as was observed for the Kodor plume. Moreover, generally distances between deltaic branches within one deltaic system are smaller than distances between estuaries of neighboring rivers. As a result, interactions between neighboring plumes formed by different rivers generally occur only during high discharge periods [86], while similar interactions between plumes formed by different deltaic branches is a permanent or almost permanent process at many World regions. However, despite a large number of deltaic rivers inflowing to the World Ocean, we are aware of only one related study that was focused on the interaction between the buoyant plumes formed by different deltaic branches of the Pearl River Delta [100].
The Kodor River has very large intra-day and synoptic variability of discharge rate due to morphology and weather conditions at its drainage basin. This variability of discharge rate induces variability of spatial extents of the Kodor plume and residence time of freshened water within the plume. As a result, the Kodor plume formed during high discharge can have different spatial and thermohaline characteristics from those formed during low discharge. In case of abrupt decrease of river discharge rate, the relatively large and mixed residual plume (formed during high discharge period) interacts with the small and freshened emergent plume (formed during the subsequent low discharge period). We report distinct frontal zones and differences in dynamics between the residual and the emergent parts of the Kodor plume. Several previous studies addressed response of river plumes to variable discharge rates [101][102][103][104][105][106][107], but limited attention was paid to interaction between parts of an individual river plume formed during different discharge conditions [108]. This feature can strongly affect spreading and mixing of freshwater discharge from small rivers in many World regions and should be considered in the related studies.
Several studies addressed interaction between coastal bathymetry and bottom-advected river plumes, which occupy the whole water column from surface to seafloor and, therefore, experience intense bottom friction [109][110][111]. In these numerical studies, river plumes were spreading over sea areas with idealized bathymetry, which was steadily sloping in the cross-shore direction and was homogenous in the alongshore direction. Influence of realistic bottom topography on surface-advected river plumes was described by Korotenko et al. [112]. Bottom-generated turbulent mixing induced by coastal circulation penetrates upward and reaches surface layer over shallow zones, therefore, increased local mixing of river plumes occurs at these zones. We presume that a similar mechanism induced intensified mixing of the Kodor plume over the shoal revealed by aerial imagery and in situ measurements. Moreover, we observed that the intense flow of the Kodor plume over this small shoal results in formation of large area within the plume with elevated salinity, which is bounded by the distinct frontal zone. We are not aware of any work describing this effect at river plumes, however, it can be typical for many small plumes with small vertical scales flowing over bathymetric features.
In this study, we address several important dynamical features of small river plumes. Aerial remote sensing revealed a quick motion of the Kodor plume border (~0.5-1 m/s) entrained by the rotating coastal eddy. Such extremely rapid response of a river plume to coastal sea circulation has not been reported before, to the extent of our knowledge. The previous studies showed that general spreading patterns of small plumes are governed by wind forcing, while the impact of ambient circulation was regarded as negligible [113][114][115]. We demonstrate that energetic features of coastal circulation, e.g., eddies, can induce high velocity motion of plume fronts and, therefore, influence dynamics of a small plume, albeit locally and during short-term periods.
The rotating eddy generated high-frequency internal waves that were propagating within the Kodor plume and dissipated at its border with the ambient sea. Aerial remote sensing also observed multiple long internal waves propagating within the Kodor plume towards the coast, as well as generation of high-frequency internal waves near the mouth of the Bzyp River and their propagation within the Bzyp plume towards the open sea. Internal waves are common features of river plumes in non-tidal seas and their surface manifestations observed by satellite imagery were reported in several previous studies [116,117]. These internal waves can significantly affect mixing of small plumes with subjacent saline sea [33]. In this study, we demonstrate the efficiency of aerial remote sensing in observations of surface manifestations of internal waves, and the ability of aerial remote sensing (in contrast to satellite observations) to measure their spatial and dynamical characteristics and to identify mechanisms of their generation.
Finally, in this study, we address the undulate structure of the sharp borders of the Kodor and Bzyp plumes that were previously observed and reported at other small plumes [52,53,118,119]. Horner-Devine and Chickadel [53] associated formation of the lobe-cleft structures observed at the Merrimack plume with subsurface vortexes that were propagating from the inner part of the plume towards its border with saline sea. Based on processing of aerial video records, we reconstruct surface circulation at the undulate fronts of the Kodor and Bzyp plume and detected similar vortexes within the lobes. However, we observed absence of vortexes outside the frontal zones, i.e., no vortexes were propagating from the inner parts of the plume towards their borders. On the opposite, we observed the recurrent process of formation, merging, and dissipation of lobes, that was not described before. Based on these results, we suggest an alternative mechanism of formation of the undulate fronts caused by baroclinic instability between the plume and the ambient sea and ballooning of local convex segments of the frontal zone in response to its small perturbations. This mechanism is in a good agreement with the reconstructed vortex circulation within the lobes and explains absence of vortexes in the inner parts of the plumes. We reveal intense transport of saline sea water across the undulate plume border as a result of merging of lobes and mechanical trapping of spots of saline sea inside the plume. It can be an important mechanism of mixing between the plume and the saline sea and should be considered together with shear-induced mixing of the plume and the subjacent sea. Satellite imagery reveals that undulate frontal zones are, therefore, the related mixing mechanism are typical for many small plumes in the World. Therefore, study of this mechanism is important in context of transformation and dissipation of freshwater discharge in the sea.

Conclusions
In this work, we focused on small buoyant plumes formed by the Kodor and Bzyp rivers located at the northeastern part of the Black Sea. We used quadcopters equipped with video camera to perform aerial remote sensing of these river plumes, which was accompanied by synchronous in situ measurements in the sea. Using an optical flow approach, we reconstructed surface velocity fields within these plumes from the obtained aerial video records. Based on aerial imagery and video records, the reconstructed surface currents, as well as in situ salinity, turbidity, and velocity measurements, we obtained new insights into spatial structure, short-temporal variability, and dynamical features of small river plumes, which are not typical for plumes formed by large rivers.
Based on the obtained aerial and in situ data, we address several different issues, including the methodology and value of the aerial observations of small river plumes, the differences between small and large plumes, the influence of multiple freshwater sources on the structure of a small plume, the influence of bathymetry features on the structure of a small plume, the interaction between small plumes and coastal circulation, the presence of internal waves in river plumes, and the presence of small-scale instabilities along the plume front boundary. The main results obtained in this study are the following. We describe strongly inhomogeneous structure of small plumes, as compared to large plumes. We suggest a new mechanism of mixing of a small plume with ambient sea as a result of baroclinic instability at its outer boundary. We describe internal waves formed within near-and far-filed parts of small plumes, which can strongly influence its mixing with ambient sea. These results are important for understanding the fate of freshwater discharge from small rivers and the related transport of suspended and dissolved river-borne constituents in many coastal sea areas in the World Ocean.
Usage of quadcopters provides ability to perform low-cost aerial remote sensing of coastal sea areas and continuously observe surface manifestations of many coastal processes. In this study, we demonstrate its efficiency in observations of small river plumes characterized by high color contrast with ambient sea, energetic motion, and high short-temporal variability. Aerial imagery can be used for visual detection and tracking of many other processes at small spatial (from meters to kilometers) and temporal (from seconds to hours) scales, which are visible neither from shipboard nor satellite imagery. Spatial scales and motion speeds of the observed processes can be reconstructed from aerial imagery and video records (Supplementary Materials). Therefore, aerial drones can provide quantitative measurements of distances and velocities at sea surface. Finally, aerial remote sensing can be very useful for operational organization of in situ measurements during field surveys, in particular, for selection of places for water sampling and hydrological measurements according to real-time position of the observed sea surface processes. As a result, future studies based on imagery and video records of ocean surface acquired from aerial drones (considering certain important limitations of their usage) and supported by in situ measurements hold promise to significantly improve understanding of various upper ocean features and dynamics.
Funding: This research was funded by the Ministry of Science and Higher Education of Russia, theme 0149-2019-0003 (collecting and processing of in situ data), the Russian Science Foundation, research project 18-17-00156 (collecting and processing of aerial imagery, study of spreading of river plumes) and the Russian Ministry of Science and Higher Education, research project 14.W03.31.0006 (study of submesoscale dynamics of river plumes).