Reach-Scale Map** of Surface Flow Velocities from Thermal Images Acquired by an Uncrewed Aircraft System along the Sacramento River, California, USA

Abstract

An innovative payload containing a sensitive mid-wave infrared camera was flown on an uncrewed aircraft system (UAS) to acquire thermal imagery along a reach of the Sacramento River, California, USA. The imagery was used as input for an ensemble particle image velocimetry (PIV) algorithm to produce near-continuous maps of surface flow velocity along a reach approximately 1 km in length. To assess the accuracy of PIV velocity estimates, in situ measurements of flow velocity were obtained with an acoustic Doppler current profiler (ADCP). ADCP measurements were collected along pre-planned cross-section lines within the area covered by the imagery. The PIV velocities showed good agreement with the depth-averaged velocity measured by the ADCP, with $R^2$ values ranging from 0.59–0.97 across eight transects. Velocity maps derived from the thermal image sequences acquired on consecutive days during a period of steady flow were compared. These maps showed consistent spatial patterns of velocity vector magnitude and orientation, indicating that the technique is repeatable and robust. PIV of thermal imagery can yield velocity estimates in situations where natural water-surface textures or tracers are either insufficient or absent in visible imagery. Future work could be directed toward defining optimal environmental conditions, as well as limitations for map** flow velocities based on thermal images acquired via UAS.

1. Introduction

There is continued interest in the development and application of remotely sensed or non-contact methods for characterizing hydrologic parameters at various scales [1]. For example, in the far-field, the recently launched Surface Water Ocean Topography (SWOT) mission is addressing spatial and temporal gaps in hydrologic monitoring using satellite radar altimetry to measure water surface elevation, slope, and extent at a global scale [2]. In the near field, ex situ sensors including radars [3] and cameras [4,5] can be deployed from both fixed and mobile platforms for detecting water level, measuring surface flow velocity, and calculating streamflow. Camera-based methods for characterizing surface velocity have seen increased adoption and widespread use by river scientists due to inexpensive hardware and proliferation of freely available image velocimetry algorithms [6,7]. Practitioners can choose from an array of algorithms including particle image velocimetry (PIV) [8,9], particle tracking velocimetry (PTV) [10], and space-time image velocimetry (STIV) [11,12]. Algorithms are increasingly being integrated into software packages designed to be user friendly and efficiently handle all phases of an image velocimetry workflow [13,14,15,16,17]. However, because cameras are passive sensors, environmental conditions can compromise the quality of imagery and ultimately the success and reproducibility of velocimetry. For example, frames used as input for image velocimetry algorithms can be susceptible to ambient lighting and shading effects that obscure measurement regions or reflected solar radiation that saturates the detector array. In some cases, a non-uniform concentration or lack of tracer particles or water-surface textures can prevent the image velocimetry algorithm from producing meaningful results.

To mitigate issues arising from insufficient tracer particles or textures, researchers have seeded the flow with particles [18,19,20,21,22] or explored alternatives to visible wavelength cameras. Infrared cameras can capture the motion of naturally occurring thermal features at the water surface. A subtle contrast in temperature at the water surface stems from turbulence that originates at the riverbed and creates boil vortices which upwell and influence the temperature of the thin (1 mm) surface layer visible to infrared cameras [23]. Although these thermal features deform as they advect with the river flow, the textural pattern can be tracked by image velocimetry algorithms. To capture these patterns, long-wave infrared cameras (8–14 μm) have been deployed from riverbanks or masts [23,24,25,26]. Cooled mid-wave infrared cameras (3–5 μm) with increased sensitivity have been used from bridges and fixed elevated platforms [27,28].

Aerial platforms provide the ability to conduct reach-scale map** of river channels. Cameras mounted on these platforms can collect near-nadir images, which minimizes the influence of perspective distortion created by viewing the river from an oblique angle. In the context of velocity map**, aerial acquisitions can provide an efficient means of collecting continuous coverage and performing synoptic investigations of river processes. For example, comprehensive spatially explicit maps of river velocity can be used to understand hydraulic controls on habitats or assist in predicting the fate and transport of contaminants. Crewed aircraft platforms have been used to collect thermal [29] and visible [30] imagery for continuous river velocimetry over large spatial extents. Recently, drones or uncrewed aircraft systems (UAS) have increased in popularity, enabling a variety of sensors to be deployed at the river-reach scale from low-cost platforms [31,32]. Although the potential of UAS-based thermal payloads for image velocimetry has been demonstrated [33,34], the authors are unaware of any continuous river-reach scale map** applications in the literature.

We designed this study to address two objectives. The first objective was to introduce an innovative sensor payload, the River Observing System (RiOS), that contains a suite of sensors including a sensitive mid-wave infrared camera for image velocimetry. The second objective was to build upon and extend the work of Kinzel and Legleiter [33], wherein UAS-based thermal image velocimetry was evaluated at two transect locations in the context of streamflow measurement. This work expands the thermal velocimetry technique to a sequence of hovering observations that enabled velocity map** over a 1 km reach in a river approximately 100 m wide. We compare the results of the thermal image-derived velocity estimates to field measurements using a conventional acoustic Doppler current profiler (ADCP). Additionally, we evaluate the spatial patterns of surface velocity captured by flying the same reach with RiOS on consecutive days. These repeat observations provide a measure of uncertainty and a means to evaluate the reproducibility of the reach-scale map** approach.

2. Materials and Methods

2.1. Study Area

The study was conducted 6–10 November 2023 along a 1 km reach of the Sacramento River located near Glenn, California, USA (Figure 1). The streamflow in the Sacramento River is influenced by upstream reservoirs, diversions, and irrigation return flows. The study reach is located approximately 46 km upstream of a U.S. Geological Survey (USGS) streamgage (#11389500 https://waterdata.usgs.gov/monitoring-location/11389500/#parameterCode=00060&period=P7D&showMedian=false, accessed on 5 May 2024. Sacramento River at Colusa [35]). The mean annual flow at the USGS Colusa streamgage, based on 77 years of record, is 317 m³/s. The study reach was used to investigate passive depth retrieval from UAS-collected multispectral imagery [36], and hydraulic models were constructed to facilitate the development and testing of image velocimetry algorithms [37]. The wetted channel width of the study reach varies between approximately 110 m at the upstream end to 150 m at the downstream end. The water-surface slope of the reach is approximately 0.0004, and the riverbed is composed of sand and gravel. Large wood present along the reach creates flow obstructions and leads to wake-induced patterns of deposition. The mean daily streamflow at the USGS Colusa streamgage decreased from 123 m³/s on 6 November to 111 m³/s on 10 November [35] (Figure 2). We continuously recorded air and water temperature with sensors that were installed along the bank within the study reach. Pre-planned cross-section lines oriented perpendicular to the primary downstream flow direction were defined at 150 m intervals. These lines provided reference locations for both the UAS-based image collection as well as the ADCP field measurements. The cross sections were designated by their stationing in meters from an arbitrary starting location upstream. For this study, we focused on eight cross sections starting at 900 m and extending downstream to 1950 m (Figure 1).

2.2. Field Measurements

To assess the accuracy of image-derived velocity estimates, field measurements of flow velocity were obtained on 6 and 7 November using a SonTek M9 (San Diego, CA, USA) ADCP enabled with a real-time kinematic global navigation satellite system (GNSS) receiver (Figure 1 and Figure 2). The manufacturer reported the accuracy of the M9 to be ±0.2 cm/s [38]. The ADCP and GNSS antenna were mounted to a pole attached to the side of a jet boat. The front of the ADCP transducer was set at depths of 0.24 m and 0.29 m below the water surface on 9 November and 10 November, respectively. ADCP measurements were collected along the cross-section lines using the RiverSurveyor Live software package [39]. At each of these cross-section locations, multiple ADCP passes were made, and the data were then processed using the Velocity Map** Toolbox (VMT) [40] to produce mean cross sections using horizontal and vertical grid node spacings of 1 m and 0.4 m, respectively. VMT does not extrapolate velocities to the water surface. The depth-averaged velocities from the mean cross sections computed in VMT are shown in Figure 1. These data are available in [41]. The deviation from the quasi-uniform spacing seen near transect 1500 is due to the presence of large wood and shallow submerged bars along the pre-planned cross-section line. These obstacles required the boat operator to adjust the ADCP collection path to avoid these hazards.

2.3. UAS-Based Image Collection

Aerial imagery was collected with a multisensor UAS payload developed by the USGS and the National Aeronautics and Space Administration (NASA) (Figure 3). The payload, called the River Observing System (RiOS), includes a thermal infrared camera, visible camera, laser range finder, inertial navigation system, embedded computer for storing data, and wireless link for transmitting data to a ground station. The infrared camera in the RiOS payload is the same camera that was used in [33], with the only difference being that the lens in the present study has a wider field of view. Specifications of the infrared camera in RiOS are listed in

Table 1. Infrared camera specifications.

Model	ICI Mirage 640 P-Series [45]
Lens	11.2 mm
Detector	Cooled Indium Antimonide
Pixel dimensions on focal array	15 μm
Number of pixels	640 × 512
Wavelength range	3–5 μm
Noise-equivalent temperature difference	<0.012 °C at 30 °C
Bits per pixel	14
Camera dimensions	111 × 96 × 131 mm
Camera weight	<765 g (without lens)
Power	12 V

. The embedded computer runs a version of the Robot Operating System (ROS) [42]. The ROS provided a convenient means of acquiring and synchronizing data from multiple sensors. Additionally, a PIV package has been developed to run in ROS, with the goal of enabling real-time image velocimetry [43]. The power supply to RiOS is contained on the payload, and all sensor systems are independent of the UAS platform. The RiOS payload enclosure is a custom fabricated carbon fiber shell that is rigidly mounted to the frame of a DJI Matrice 600 Pro hexacopter (Shenzhen, China) [44] using a “toad in the hole” quick-release mechanism. Although incorporating a gimbal mount to improve image stabilization was considered, the maximum weight limitation of the platform prevented this design [44].

A flat and level area in a field adjacent to the river was selected for UAS takeoff and landing. Map** flights were conducted between 9:00 and 10:30 a.m. local time on 9 November and between 8:00 and 8:30 a.m. local time on 10 November (Figure 2). The winds were generally calm in the mornings when the flights were conducted. For each flight, the pilot used a tablet running the Ignis application [46] to navigate to waypoints collocated with the pre-planned cross sections using the internal GPS of the UAS. At each of these waypoints, the pilot flew the UAS to a height of approximately 250 m above the center of the river channel, an altitude from which the field of view was sufficient to capture both riverbanks within the images. A wireless data link between the UAS and ground station facilitated near real-time viewing of the thermal imagery by the pilot. Using this data link, the pilot was able to adjust the position of the UAS to center the river channel between the banks as captured in the image frame. This level of spatial coverage ensured that adequate features along the banks were included in the thermal imagery. These features were used in the image stabilization and geo-registration steps of the processing workflow discussed in Section 2.4 and shown in Figure 4.

Once positioned at the waypoint location and with both banks present in the live image feed, the pilot hovered while maintaining the position for approximately one minute. We refer to the period spent hovering at the waypoint location as the dwell. The dwell was based upon a previous examination of the effect of dwell on ensemble PIV accuracy, which found that, provided an appropriate frame rate was selected, longer dwell times did not provide more accurate velocity estimates [47]. The thermal images were saved to the memory of the embedded computer on the payload as ROS *.bag files. These *.bag files are available in [41]. The flight time of the UAS with the RiOS payload was approximately 13 min. During this period, the UAS was able to transition from the takeoff/landing area, acquire imagery at several waypoints, and have sufficient remaining power to safely return, land, and power down so that the batteries could be replaced.

2.4. Thermal Image Processing

After each flight, the *.bag file was downloaded from the embedded computer and saved to external storage media for further processing and analysis. Foxglove Studio, a freely available software for viewing ROS *.bag files, was used to open, play, and review the files [48]. The start and end times of each dwell were identified by reviewing the imagery and noting the appropriate timestamps in the *.bag file. For each dwell in the *.bag file, MATLAB [49] was used to execute the processing workflow shown in Figure 4. The mean frame interval was calculated from the range of timestamps for all of the image frames included in the sequence. The ratio of the mean frame rate to the target frame interval, in this case 1 Hz, was calculated to determine the image skip factor necessary to achieve the specified target frame interval. The frames were then extracted from the *.bag file for each dwell using that skip factor. Finally, a target duration of 60 s was used to truncate the set of extracted frames. Ultimately, a set of 69 frames was retained for each dwell with a frame interval of approximately 0.87 s per frame. The pixel values in each frame were normalized by the maximum value within the scene, providing a scaled image with values between 0 and 1 (Figure 5).

A series of image-processing steps were applied to each frame in a given dwell using MATLAB (Figure 4). First, an interactive contrast enhancement tool available in MATLAB allowed the analyst to set the upper and lower limits of the image histogram. This step provided a way to hone in on a range of values that accentuated the features in the image (i.e., trees along the banks and large wood in the channel) necessary for the subsequent stabilization step. Stabilization removes the platform motion and ensures that the movement of thermal features can be attributed solely to the river flow. The stabilization algorithm in the workflow finds distinct features in all the frames in the dwell, identifies the correspondences between these features and those in a reference frame (typically the middle frame in the dwell), and then performs an affine transformation to align all images in the dwell sequence with the reference frame. The stabilization algorithm was robust and did not use any features within the water, so masking of the water surface was not required for this step. The next step in the workflow is background subtraction to mitigate lens or sensor artifacts that are present in every frame in the dwell, such as vertical banding. A fast Fourier transform (FFT) bandpass filter was then applied, followed by a histogram equalization. These steps enhanced the contrast of the thermal features in the river channel. Another algorithm then identified the common area of overlap between the frames in the dwell and the final region of interest defined by interactively digitizing the boundary of the active river channel. The thermal imagery was georeferenced to a base color orthophotograph that was collected during the same week (Figure 1). Corresponding features (large wood in the channel and distinct vegetation patterns along the banks) were selected from both the thermal image and the base, and an affine transformation was applied to warp the thermal frames to the base image (Figure 5).

2.5. PIV Algorithm

The processed thermal imagery had a pixel size of 0.34 m and was used as input for a PIV algorithm. We used source code from the Toolbox for River Image Velocimetry using Images from Aircraft (TRiVIA) [16] for the PIV analysis (Figure 4). TRiVIA calls the core ensemble PIV algorithm that is available in PIVlab [13]. Instead of specifying the interrogation area and step size in pixels as is performed in PIVlab, TRiVIA allows the user to define the output spacing between vectors in units of meters. An output spacing of 5 m (equivalent to approximately 15 pixels) was selected for the imagery collected along the study reach. This resolution provided approximately 20 measurements of velocity across the channel, sufficient to capture cross-stream variation in the flow and the influence of various obstacles, such as large wood, on the flow field. A single pass of the PIV algorithm provided reasonable output, and vector smoothing and infilling were applied.

2.6. PIV versus ADCP Accuracy Assessment

The ADCP field measurements of depth-averaged velocity were used in a PIV accuracy assessment algorithm that was also adapted from the TRiVIA software, Version 2.1.1 [16]. For each PIV velocity estimate, this algorithm searches for all ADCP measurements within a search radius that is set to one half the output vector spacing. In situations where the search radius included more than one ADCP measurement, all the ADCP measurements within the search area were averaged to obtain a single value for comparison to the velocity estimated via PIV. After each PIV estimate is paired with a velocity based on the field measurements, several metrics are computed: observed vs. predicted coefficient of determination ($R^2$), the slope and intercept of the regression line, normalized root mean square error (RMSE), and normalized bias. The $R^2$ value quantifies the level of agreement between all the paired values. The slope and intercept of the regression line also provide some additional information on the relationship between the ADCP depth-averaged velocity and corresponding PIV-derived surface velocity. The RMSE evaluates the differences between the observed (i.e., ADCP) and the predicted (i.e., PIV) velocity magnitudes and provides an index of precision. In calculating the mean bias, the errors are defined as PIV-ADCP. For example, a positive value of the mean bias indicates that PIV overestimated the ADCP measurements on average, and a negative value for the bias indicates that PIV underestimated the ADCP measurements on average. The bias and RMSE were normalized by dividing the mean bias and RMSE by the mean of the ADCP velocity magnitudes measured in the field.

2.7. Day-to-Day PIV Accuracy Assessment

We evaluated the agreement between the velocity fields from 9 and 10 November using the same procedure implemented in the TRiVIA software and discussed in Section 2.6 to compare PIV-derived velocities to ADCP field measurements. However, in this case, instead of using the ADCP field data, we compared the PIV-derived velocities from one day to the next. In addition, for these day-to-day comparisons, we computed several metrics that considered both the velocity magnitude and the vector orientation using concepts introduced in [50]. The first of these metrics, the weighted relevance index (WRI), quantifies how well two vector fields are aligned. It uses a scheme that assigns more weight to misaligned vectors with higher velocities. The median velocity magnitude of each vector field is used for normalization, which ensures an equivalent treatment of each vector field. For example, if the vector fields from 9 and 10 November were perfectly aligned, the WRI would be 0, whereas a WRI of 0.5 would indicate a 90° misalignment between the two vector fields. The weighted magnitude index (WMI) provides an index to assess the agreement between two velocity fields based only in terms of magnitude. It is computed by dividing the difference in vector magnitude at a given location by a median value that includes every magnitude in the two velocity fields. If the vectors for 9 and 10 November have the same magnitude, the WMI will be 0. If there are large differences in the magnitudes, the WMI will be >1. Finally, because good agreement in terms of magnitude does not necessarily ensure good agreement in terms of orientation and vice versa, a combined summary metric, the combined magnitude and relevance index (CMRI), is computed as the average of the WRI and WMI.

3. Results

The PIV versus ADCP accuracy assessment procedure outlined above yielded results indicating good agreement between the ADCP field measurements and PIV estimates; $R^2$ values ranged from 0.59 to 0.96 on 9 November and from 0.60 to 0.97 on 10 November (

Table 2. Comparison of PIV and ADCP velocities along transects for 9 November and 10 November. Transect (XS), Observed (Obs.), Predicted (Pred.), Normalized (Norm.).

XS	Obs. vs Pred.	Norm.	Norm.	Slope	Intercept	n
	${\mathit{R}}^{\mathbf{2}}$	RMSE	Bias
	9 November 2023
900	0.94	0.12	0.10	1.16	−0.04	14
1050	0.95	0.05	0.02	0.96	0.04	19
1200	0.87	0.15	0.11	1.03	0.05	17
1350	0.96	0.06	0.04	1.05	−0.01	22
1500	0.86	0.14	0.11	1.28	−0.16	25
1650	0.59	0.20	0.18	1.37	−0.20	23
1800	0.68	0.19	0.16	1.27	−0.11	22
1950	0.68	0.23	0.01	2.25	−1.34	17
	10 November 2023
900	0.75	0.14	0.04	1.02	0.01	16
1050	0.97	0.04	0.00	1.00	0.00	18
1200	0.94	0.11	0.04	1.14	−0.06	20
1350	0.72	0.09	−0.02	0.84	0.10	20
1500	0.77	0.17	0.14	1.18	−0.04	22
1650	0.83	0.18	0.18	1.43	−0.25	23
1800	0.60	0.17	0.14	1.04	0.11	21
1950	0.86	0.09	0.06	1.36	−0.33	15

). The normalized RMSE ranged from 0.05 to 0.23 on 9 November and from 0.04 to 0.18 on 10 November, implying that on average, the PIV estimate was within 4 to 23% of the mean depth-averaged velocity measured in the field. The normalized mean bias was generally positive for both days, which indicated that overall, the PIV surface velocity magnitudes were greater than the corresponding ADCP depth-averaged velocity values. This result was expected because theoretical velocity profiles (e.g., power law, logarithmic) indicate that the surface velocity is greater than the depth-averaged velocity [51]. The regression line slopes provided some additional information on the relationship between the surface velocity and depth-averaged velocity at each of the transects in the reach. We also aggregated the transect-level comparisons to create a reach-level comparison for each day; these plots are shown in (Figure 6).

Ideally, the above comparisons would be made for contemporaneous ADCP measurements and PIV estimates. Since the ADCP velocity measurements were made several days prior to the PIV estimates, at a flow rate that was 10% greater, some uncertainty is associated with these comparisons. In any case, although a detailed examination or fitting of ADCP velocity profiles is beyond the scope of this manuscript and has been treated elsewhere [52], the slopes of the reach-level comparisons in Figure 6 can be used to obtain estimates of the velocity index $\alpha$, which is defined as the ratio of the depth-averaged velocity to the surface velocity. The inverse of the regression line slope provides an approximation of $\alpha$ and was 0.862 on 9 November and 0.854 on 10 November. These values are typical of those found in the literature and close to the widely used default $\alpha$ of 0.85 [53]. Kinzel and Legleiter [33] showed that by assuming $\alpha =0.85$, thermal PIV and ADCP depth-averaged velocities compared favorably in a previous study on the Blue River, Colorado, USA.

Similar to the decreasing streamflow measured at the USGS Colusa streamgage, the water temperature steadily declined over the course of the week. Figure 2 shows the air and water temperatures at the time of the UAS image collections at the study site. The air temperature exhibited a typical diurnal cycle. The thermal imagery collected on 9 November was collected later in the morning than the imagery on 10 November. During the later time on 9 November, the air temperature was as much as 6 °C warmer than on the morning of 10 November. Although in some cases the imagery from 10 November seemed qualitatively better, in terms of improved sharpness of the thermal contrast within the river, there was not a clear difference between the PIV versus ADCP velocity metrics from one day to the next on either a transect-by-transect basis or at the reach level. Instead, regardless of the day the survey was conducted, variations in the metrics were observed from one cross section to the next. These variations were more likely a consequence of the different ranges of velocities sampled over the cross sections and not due to the image quality alone. The narrower upstream portion of the reach (cross sections 900 through 1350) had relatively slower ADCP-measured velocities compared to the downstream, wider portion of the reach (cross sections 1500 to 1950) (Figure 1). These higher velocities in the downstream section could be attributed to the presence of large wood obstructions and emergent bars that concentrated the flow and also led to shallower depths. Given these different velocity regimes and the several-day interval between the field measurements and UAS flights, the PIV-derived velocities showed strong agreement with the ADCP velocities (Table 2 and Figure 6).

The surface velocity fields computed via PIV for each hover observation were combined to produce reach-scale velocity maps for both 9 and 10 November (Figure 7). Although the surface velocity of the study reach is characterized with high spatial resolution, there are gaps between the vector fields derived from the various dwells. The gaps are due to the altitude at which the imagery was acquired and the corresponding field of view of the camera, as well as the spacing of the cross-section observations that were selected. The gaps can also be attributed to the crop** we performed to avoid edge effects on the upstream and downstream margins of each image sequence. The upstream four rows of PIV-derived velocity estimates and downstream three rows were omitted from the reach-scale plots. We had to apply a buffer at the upstream edge because the features being tracked by the PIV algorithm could not be tracked by the algorithm until they developed and entered the scene. In spite of the gaps in coverage, the plots show a similar spatial pattern to the ADCP measurements, with slower velocities measured in the upstream, narrower portion of the reach (cross sections 900 through 1350) followed by higher velocities in the downstream portion (cross sections 1500 to 1950). Both plots show a high velocity thread develo** at cross section 1500 in the center of the channel, and these higher velocities are maintained through cross section 1650. Along cross section 1800, an emergent bar splits the flow, with a high velocity core to the south (river left). In the final cross section, 1950, the same high velocity core is directed along the south half of the channel with an area of much more stagnant flow present along the river right.

A comparison of the velocity magnitudes from 9 November and 10 November using the accuracy assessment procedure outlined in Section 2.6 produced the following results: an $R^2$ of 0.86, a normalized RMSE of 0.17, and a normalized bias of −0.02. Figure 8 shows the relationship between the velocity fields collected on 9 November and 10 November. Across most of the range of observed velocities, there is good agreement. However, more scatter is seen at the lower end with slower velocities from 10 November being paired with faster velocities on 9 November. This could be potentially attributed to the georeferencing process, which could have introduced mismatches along the banks and near the obstacles to the flow.

The WRI, WMI, and CMRI calculated by comparing the velocity fields from one day to the next are provided in Figure 9. The WRI shown in the top panel of Figure 9 indicates that the velocity fields from 9 November and 10 November were well aligned with the histogram showing that most values were less than 0.01. The map on the top right shows the spatial distribution of WRI. Except for one data point of locally high WRI in hover 1950, the WRI is generally near zero throughout the entire domain. The WMI shows more variability, with the highest values of WMI along the channel banks where some variation in PIV results would be expected because of edge effects. Additionally, some locally high values of WMI are seen near some of the exposed bars and large wood. This is to be expected because the flow field is more turbulent in these areas, resulting in a greater range of velocity magnitudes. The CMRI highlights many of the same areas as the WMI and for the same reasons described above.

4. Discussion

Image velocimetry in rivers can be challenging when natural tracers or surface textures are absent, sparse, or cannot be ubiquitously resolved using conventional visible cameras. Distributing natural or artificial seeding material poses logistical difficulties for large rivers, and the effects of seeding density and clustering on image-derived velocity estimates require careful consideration [19]. These impediments to image velocimetry were particularly acute in the Sacramento River at the time of our field study. The streamflow of the Sacramento River was approximately 115 m³/s, which is substantially lower than the mean annual flow of 317 m³/s. With this relatively low flow rate, natural floating material was lacking. Although water surface textures downstream of flow obstructions could be observed obliquely from the ground, imaging these features with the visible camera on the RiOS payload was difficult because of the nadir viewing geometry and low sun angle at the time of our UAS collections. However, the images collected by the RiOS thermal camera provided an alternative means for image velocimetry. The thermal contrast at the water surface could be reliably tracked, and the velocity fields from individual hover acquisitions could be combined into reach-scale maps.

This manner of comprehensive and continuous map** would not have been possible with in situ sensors such as ADCPs. Velocity map** with ADCPs can be time consuming and disruptive, with the potential of leaving large gaps in portions of the channel which cannot be safely or efficiently accessed. Submerged and emergent obstacles in particular can pose a significant hazard risking injury to personnel, damaging submerged acoustic equipment mounted to the boat, or even the boat itself. For example, areas of the Sacramento River had to be avoided during our ADCP collections because of the risk of grounding the boat on shallow areas or striking submerged or emergent large wood.

At present, thermal imagery introduces some additional steps beyond those included in a conventional image velocimetry workflow. Although visible imagery generally requires some pre-processing operations such as contrast limited adaptive histogram equalization (CLAHE) [54] to accentuate image features prior to PIV, we found some manual adjustment of our thermal imagery was necessary to accentuate textures. Eltner et al. [31] also used image processing to enhance the signal-to-noise ratio in their thermal imagery. The analyst performed contrast enhancement following collection, which was a prerequisite for the stabilization step of our workflow. This need for manual intervention has implications for automating the image pre-processing step onboard the UAS [43] and computing PIV estimates in real-time. Conducting further testing and comparison of the results of automated algorithms with and without user intervention would be informative.

A past deployment of this thermal camera along the Blue River in Colorado, USA, on a UAS [33] was conducted at 06:50 a.m. local time, when the air temperature was cold (−5.3 °C) and the water in the river was relatively warm (6.9 °C). The thermal imagery under these conditions showed good contrast at a flying height of 100 m, and PIV was successful following image enhancement. In a related thermal image velocimetry study, Schweitzer and Cowen [28] presented hourly images collected from 18:00 to 09:00 using a mid-wave infrared camera with higher resolution than the camera in the RiOS payload. Their camera was elevated with a forklift above a tributary to the Sacramento River. They noticed decreased sharpness of their imagery at both the beginning and end of this period. Like our Blue River study [33], the water temperature in the tributary was higher than the corresponding air temperature. The authors noted that when the difference between the wet bulb temperature (which includes both air temperature and humidity) and the water temperature was greater than 3 °C, image quality improved. For our present study on the Sacramento River, to achieve a condition where the air temperature was colder than the river water would have required setting up and conducting our flights in low-light conditions. Similarly, scheduling our collections at night to achieve this temperature differential would have added complexity and raised regulatory and safety concerns. Realistically, neither the early morning nor overnight periods are practical for conducting operational UAS measurements. In our present study, the air temperature was greater than the water temperature in the Sacramento River for both 9 November and 10 November. Although, qualitatively, the imagery was sharper on 10 November than 9 November, after image enhancement, PIV was successful for both data sets, and the results did not differ markedly from one day to the next.

Although an absolute temperature difference between the air and the water is required to produce thermal contrast, both turbulence within the water column and the dynamic response of water temperature to solar radiation over the course of a day are likely to affect the level of thermal contrast. The influence of humidity on image quality identified by Schweitzer and Cowen [28] could also be explored further. Experimenting with the RiOS thermal camera over a broader range of conditions could improve our understanding of the minimal and potentially optimal operational conditions for thermal PIV and help outline practical limitations for use.

There are several planned enhancements to RiOS that would create additional capabilities and could open up additional directions for research. Our present concept of operations involves relying on reference points or features on the ground to provide both a means to stabilize frames and also for georeferencing. Although georeferencing places the imagery in the correct spatial position and is necessary for comparison to GNSS-enabled field measurements and reach-scale map**, the georeferencing step in and of itself is not necessary for performing PIV. Our present field protocol involves a camera calibration step. Prior to attachment to the UAS, RiOS is placed on its side, and a calibration checkerboard is moved through the field of view of the cameras. The checkerboard pattern is made of a material that can be resolved in both the visible and thermal imagery. This calibration procedure allows the precise focal length of the cameras as well as other intrinsic parameters to be calculated. With the focal lengths of the cameras and the distance from the cameras to the water surface measured by the onboard laser rangefinder, the pixel size of the imagery can be determined. With this scaling factor and a knowledge of the camera frame rates, PIV estimates can be directly computed. However, the problem of image stabilization remains. One solution can be the transition to a UAS that can support a gimbal-mounted payload. Another potential solution would be to perform direct georeferencing of the imagery. This process would require camera intrinsic parameters, precise positioning from a real-time kinematic GNSS, and orientation data from an inertial measurement unit (IMU) mounted on the camera. With the implementation of either solution, the river banks would not be required to be in the field of view. Thus, the UAS could fly lower over the river channel and collect finer-scale velocity data. In this manner, the RiOS payload could perform real-time PIV and follow pre-planned or potentially autonomous paths to seek out the most interesting or variable regions of the flow field.

Although we used a more traditional hovering acquisition in this study, alternative flight plans could be explored. We have shown that image velocimetry is possible using fixed wing aircraft that follow a river course, rather than hovering in place [30]. Moving aircraft river velocimetry (MARV) could be attempted with a UAS programmed to follow a series of pre-planned flight lines. Provided sufficient forward overlap is maintained with short, consistent image capture intervals, imagery could be obtained more efficiently and with fewer gaps in coverage. Additionally, if radiometric calibration is applied to the thermal images, reach-scale absolute temperature maps could be produced from such a workflow. These data could be used to characterize thermal regimes for ecological processes [55] or identifying and map** groundwater–surface water exchange [56].

5. Conclusions

In this paper, we demonstrated the potential for spatially explicit velocity map** at the river-reach scale using an innovative UAS payload, RiOS, containing a sensitive mid-wave infrared camera. Small temperature differences at the water surface caused by turbulent mixing processes were captured in the thermal imagery and provided natural flow tracers in the Sacramento River. Image sequences collected while hovering over the channel were used as input to an ensemble PIV algorithm to estimate surface velocities. The image-derived velocities from multiple hovering observations were then combined to construct nearly continuous reach-scale velocity maps. The velocity estimates showed good agreement with field measurements of velocity acquired with an ADCP ($R^2$ from 0.59–0.97). Additionally, velocity maps collected on consecutive days during a period of relatively steady flow showed similar patterns. We highlight several considerations and caveats with this approach. With respect to image pre-processing, the thermal imagery required a series of steps beyond what has been used for visible images to accentuate the contrast of thermal features. At present, the payload does not use a gimbal and instead relies upon feature-based image stabilization, which dictates flying high enough to ensure that the images include both river banks. Images were tied to ground coordinates using tie points from a contemporaneous orthophotograph, which introduced a laborious and potentially subjective step in the processing workflow. Although this work demonstrated that image velocimetry is possible when the air temperature is greater than the water temperature, further work is needed to assess the quality of thermal imagery collected by the payload at both diurnal and seasonal time scales to identify the environmental conditions conducive to a favorable signal-to-noise ratio.