Detailed understanding of submesoscale processes and their role in global ocean circulation is constrained, in part, by the lack of global observational datasets of sufficiently high resolution. Here, the potential of thermosalinograph (TSG) and Visible Infrared Imager Radiometer Suite (VIIRS) data is evaluated, to characterize the submesoscale structure of the near-surface temperature fields in the Gulf Stream and Sargasso Sea. In addition to spectral density, the structure function is considered, a statistical measure less susceptible to data gaps, which are common in the satellite-derived fields. The structure function is found to be an unreliable estimator, especially for steep spectral slopes, nominally between 2 and 3, typical of the Gulf Stream and Sargasso regions. A quality-control threshold is developed based on the number and size of gaps to ensure reliable spectral density estimates. Analysis of the impact of gaps in the VIIRS data on the spectra shows that both the number of missing values and the size of gaps affect the results, and that the steeper the spectral slope the more significant the impact. Furthermore, the TSG, with a nominal resolution of 75 m, captures the spectral characteristics of the fields in both regions down to scales substantially smaller than 1 km, while the VIIRS fields, with a nominal resolution of 750 m, reproduce the spectra well down to scales of about 20 km in the Sargasso Sea and 5 km in the Gulf Stream. The scales at which the VIIRS and TSG spectra diverge are thought to be determined by sensor and retrieval noise.
All numerical models of physical processes make assumptions about physics at scales smaller than that of the model-grid spacing. This is referred to as subgrid-scale parameterization. Recent observational and numerical work show the physics in the upper ocean at scales smaller than the first baroclinic mode [generally referred to as the submesoscale and not resolved in global ocean general circulation models (OGCMs)] to be rich in structure and processes. Moreover, these processes may play a significant role in, and interact with, the larger-scale circulation, for example, through their impact on air–sea fluxes, vertical motion, and particle redistribution (e.g., Lapeyre and Klein 2006; Capet et al. 2008; Thomas et al. 2008; Molemaker et al. 2010; D’Asaro et al. 2011; Fox-Kemper et al. 2011).
To date, most studies of submesoscale processes have been numerical owing to the paucity of observational datasets at these scales. A number of multiplatform experiments have been undertaken at these scales [e.g., Scalable Lateral Mixing and Coherent Turbulence (LatMix)1; Ocean Surface Mixing, Ocean Submesoscale Interaction Study (OSMOSIS)2], but these experiments are expensive and very limited in space and time. On the other hand, until recently, routine sampling of the upper ocean, undertaken either over large spatial scales and/or large temporal scales, has been relatively coarse: the highest-resolution satelliteborne observations of the ocean surface have sampled at approximately 1 km, but the noise associated with these samples rendered them of little use at scales smaller than (Deschamps et al. 1981, 1984). Gliders have been employed to sample along repeat tracks (e.g., Todd et al. 2011; Cole and Rudnick 2012); however, these sections are still relatively isolated and coverage is far from global. Furthermore, with the exception of a small number of datasets, the spatial scale at which long-term, repeated ship surveys sample tends to be and coarser. The latter is determined by the temporal sampling rate and ship speed. Since these observations tend to be taken from commercial vessels, there is little control over the ship speed and the sampling rate is constrained by the instrument, ship-to-shore transmission rates, and/or onboard storage capacity. Sections from research vessels, repeated as the ships go to and from their home port, may be of interest for the study of coastal processes, but the repeated portion of these sections tends to be short and there is little regularity in their frequency.
The objective of the work presented herein is to determine the spatial fidelity of shipborne thermosalinographs (TSGs) and the Visible Infrared Imaging Radiometer Suite (VIIRS) carried on the recently launched Suomi–National Polar-Orbiting Partnership (Suomi-NPP) satellite. The TSGs have been deployed on a number of ships of opportunity but not widely used because of concerns with data quality. The interest in VIIRS stems from the fact that instrument noise is low, making for higher-quality sea surface temperature (SST) retrievals compared with those of other operational satelliteborne infrared instruments with similar resolution and coverage (Schueler et al. 2002). VIIRS samples the earth’s surface twice daily with a spatial resolution of 750 m for the operational SST product, compared with for Advanced Very High Resolution Radiometer (AVHRR) and Moderate Resolution Imaging Spectroradiometer (MODIS) products. Unfortunately, infrared radiometers are susceptible to cloud contamination so retrieved SST fields tend to be spatially patchy. A second objective of this project, therefore, is to explore the role of data gaps and their characteristics to better understand differences between in situ and satellite-derived SST datasets.
Different methodologies for deriving the statistical measures are compared for both observational and synthetic datasets. Specifically, we compare temperature spectra and structure functions obtained from VIIRS with those obtained from an acoustic Doppler current profiler (ADCP) and a TSG, both mounted on the Motor Vessel (MV) Oleander. Temperatures acquired as part of the ADCP system, although coarser in spatial resolution (3–5 km) than the TSG data (~75 m), are included here as a reference dataset in that their spectra have been published in recent studies of submesoscale processes in the region (e.g., Wang et al. 2010; Callies and Ferrari 2013). The Oleander route transects several distinct dynamical regimes; our analysis focuses on the Sargasso Sea and the Gulf Stream (Fig. 1). To test the reliability of the different methods in each of these regions, we also apply them to synthetic datasets with known properties.
The datasets used are introduced in detail in the next section. In section 3 we provide an overview of the different spectral methods used. Analysis of the TSG and VIIRS data is provided in section 4. We close with a summary and discussion in section 5. In three appendixes we discuss a simple model for the ADCP sensor response and apply the methods used to synthetic data with known properties and provide a list of acronyms.
a. Oleander region and subregions
The area of our study is the region between Hamilton, Bermuda, and Port Elizabeth, New Jersey, along the route of the Oleander, a container ship that makes weekly round trips between these two destinations (see Fig. 1). The ship traverses several distinct dynamical regimes: shelf, slope sea, the Gulf Stream, and the Sargasso Sea. While delineation of these regions is somewhat ambiguous because the Gulf Stream meanders significantly in the vicinity of the Oleander track, we define the Gulf Stream region as that between 36° and 39°N, and the Sargasso Sea region as that between 32.5° and 36°N. These choices ensure that the energetic Gulf Stream front is always contained in the Gulf Stream region and, conversely, that the Sargasso Sea region is relatively quiet.
b. Oleander data
Three of our datasets are derived from measurements taken by instruments mounted on the Oleander.
1) Acoustic Doppler current profiler (ADCP)
At the outset of the Oleander Project in 1992, the Oleander was equipped with a narrowband 150-kHz ADCP mounted in a sea chest in the keel (Flagg et al. 1998). The instrument was replaced by an Ocean Surveyor 75-kHz ADCP in 2004. ADCPs primarily measure velocity; however, to correct for the effect of temperature in the speed of sound, a thermistor measures temperature in the sea chest at about 5–6-m depth, depending on the load of the Oleander. In our study, we use the processed time series from the 75-kHz instrument from 1 January 2005 to 21 November 2013, which is bin averaged and has a temporal resolution of 3 min, corresponding to about 1.5 km at a typical cruising speed of 16 kt (8 m s−1).
The accuracy of the thermistor in the 75-kHz ADCP is 0.1 K and the resolution is 0.027 K. Because the sensor is glued to insulating material next to the array stack, it is expected to respond relatively slowly, but its response time is not well known.3 We estimate the response time of the ADCP thermistor to be between 5 and 10 min as described in appendix A. A slow response time such as this has implications with regard to the estimated spatial spectra as discussed in section 4a.
The Oleander ADCP data were downloaded from the Oleander data portal,4 and inspected for quality following Wang et al. (2010). A Barnes filter (Barnes 1964) with a decay scale of 1.5 km, which is similar to the resolution of the processed time series, was then used to interpolate sections onto an along-track grid with a resolution of 1.25 km to allow for application of standard algorithms, such as the discrete Fourier transform (DFT). For computational efficiency, the Barnes filter makes use of data points within a circle of radius 3 km centered on each grid point; if no valid data are available for a given grid point, then it is flagged as missing.
2) Thermosalinograph (TSG) and remote temperature sensor (TEX)
The Oleander is equipped with a thermosalinograph that measures temperature and salinity at the engine intake. The system was operated by the National Oceanic and Atmospheric Administration (NOAA) from 2001 to fall 2013 and by the Oleander Project from September 2014 to present. It consists of a Sea-Bird Electronics (SBE) 45, which measures temperature and salinity from the seawater intake in the interior of the ship. The accuracy of the temperature sensor is 0.002 K, the resolution is 10−4 K, and the response time is 0.5 s. Because temperature can be altered while water is flowing through the intake pipes (~4.5 m before the SBE 45 sample is taken), the thermosalinograph was equipped with an additional SBE 38 remote temperature sensor (TEX; accuracy: 0.0001 K, resolution: 0.00025 K, and response time: 0.5 s) from September 2007 to fall 2013. The SBE 38 measures temperature directly at the intake, that is, the temperature “external” to the Oleander. In our study we analyze the temperature data from September 2007 to June 2013 from the SBE 45 (TSG) and SBE 38 (TEX) sensors. Analyzing both time series provides the opportunity to explore whether the quality of the “interior” measurements is sufficient to provide reliable estimates of spectral slopes. The data were downloaded from the NOAA Atlantic Oceanographic and Meteorological Laboratory (AOML) website, where links to the datasets were made available upon request. The sampling frequency of the TSG and TEX sensors is one measurement per 10 s, such that their horizontal resolution (~75 m at a speed of 16 kt) is significantly finer than that of the ADCP.
For quality control, we use the standard quality control by NOAA with regard to whether a grid point is located over water or land. In addition, we perform our own quality control (QC) as part of the project. Specifically, we test whether temperature is in a reasonable range (0°–33°C). We also eliminate spikes from the time series; spikes being defined as data points with temperatures at least 0.5 K above (or below) those from both adjacent data points. A particular problem with the TSG data is that sometimes the pump does not operate properly. Data measured at these times are characterized by higher-than-normal temperatures and very low variability. We eliminate such faulty data through manual inspection. As with the ADCP data, the TSG and TEX sections are interpolated on an equidistant along-track grid with a resolution of 75 m, using the length scales m and m.
c. Visible-Infrared Imager-Radiometer Suite (VIIRS)
Our final dataset consists of SST fields derived from the VIIRS instrument launched on the Suomi-NPP satellite in October 2011. The SST product used here [Joint Polar Satellite System (JPSS) VIIRS sea surface temperature environmental data record (EDR)5] was retrieved from the VIIRS “moderate resolution bands,” which have a nadir resolution of 750 m and, because of the way in which the instrument is configured, decreases very slowly to approximately 1600 m at the scan edge, a ground distance of approximately 1500 km from nadir (Seaman et al. 2014; Schueler et al. 2013). For this study, we use only the best-quality data, quality level 1. This removes most clouds as well as a number of pixels that were cloud free. To compare the VIIRS data to the Oleander sections, we use an isotropic Barnes filter with m to interpolate SST on an equidistant grid with a resolution of 750 m along a nominal Oleander section with start points and endpoints at 64.7°E, 32.45°N and 72.85°E, 39.975°N, respectively, and a constant latitude–longitude angle. For each grid point, data points within a distance of (1.5 km in this case) of the point of interest are considered. When no valid data are available within , the data point is flagged as missing.
Different methodologies are used throughout the literature to estimate spectral slopes. Here, we focus on two methods: The most widely used is the DFT, generally determined by a fast Fourier transform (FFT) algorithm designed to efficiently determine the DFT (e.g., Wang et al. 2010). The second approach is based on the second-order structure function (Deschamps et al. 1981; Wald 1983; McCaffrey et al. 2015). Here we review the underlying theory and describe the steps involved in each of the methods, and introduce two measures for rating the quality of data with missing points.
a. Relation between spectral density, autocorrelation function, and structure function
In the case of an infinite, continuous time series , the autocorrelation function is the Fourier transform of the spectral density and is linearly related to the second-order structure function. Consequently, when the spectral density follows a power law , the structure function has the from , and the two slopes are related by
(e.g., Deschamps et al. 1981; Frisch 1995). In the case of discrete finite time series with , as is the case in the present study, the relations between P, , and differ. The autocorrelation function is defined by
with . The discrete Fourier transform of is
where , with being the total length of the record, and is its temporal resolution. The power spectral density is given by
where denotes the complex conjugate of . To illustrate the relation between spectral density and autocorrelation, we substitute for index m, which gives
Substitution of (2) and using the symmetries of the exponential function gives
where is the biased autocorrelation function. The second-order structure function is related to the autocorrelation function by
where the expected value of the last term in (7) is for a homogenous time series.
Equation (6) illustrates that in the case of a discrete finite time series, the power spectral density is related to the biased autocorrelation function, rather than , as would be the case for an infinite, continuous time series . Moreover, dynamical regimes with constant slopes only extend over a finite range of scales in the ocean. As a result, the precise circumstances under which (1) adequately describes the relation between slopes of the structure functions and spectra are still under debate (e.g., Webb 1964; Huang et al. 2010; McCaffrey et al. 2015). In this study we test the accuracy of relation (1) for the VIIRS and TSG datasets.
All our datasets have gaps. Because the VIIRS measurements are strongly affected by cloud cover, they have significantly more missing data along the Oleander track than the in situ sections. As a result, it is impossible to select a sufficient number of (nearly) complete VIIRS sections for our analysis. This is potentially problematic because data gaps, depending on their nature (i.e., their number and length), affect the estimate of spectral densities. Thus, for selecting the VIIRS sections used, we seek an approach that balances the number of available sections with the possible error due to the gaps contained therein. Here, we introduce two measures to characterize the gaps. These measures are used in our study to assess the quality of a section and its potential for improving the ensemble mean.
Consider a time series with N data points , with being the data points with valid values, and being the data points without valid values. Term is the number of valid points, and is the number of points with missing data. We define
as the ratio of valid to total data points in each section, and we define the “noise cohesion” as
where is the number of neighboring data-point pairs , where neither nor is a valid data value; that is, both belong to (e.g., Cayula and Cornillon 1992). The second equality in (9) illustrates the (nonlinear) relationship between the noise cohesion and the average length of the gaps, , where is the number of gaps. So, for any two of N, Q, and fixed, C increases toward 1 as the remaining parameter increases (note that and ). Our results (discussed below) suggest that C provides a useful measure for obtaining a first-order characterization of the gap structure. Figure 2 shows 2D histograms, indicating the relative abundance of VIIRS sections as a function of Q and C in the Sargasso Sea and the Gulf Stream.
1) Discrete Fourier Transform (DFT)
A common method to estimate spectral densities is to use FFT algorithms to obtain the Fourier transform and from this the spectral density. These algorithms require gap-free time series6 on equidistant grids; that is, gaps, if they exist in the original time series, must be filled prior to applying the algorithm. Gaps result from cloud cover in the satellite-derived SST fields, and from bubbles or intermittent system failures in the in situ time series. Here, we use a Barnes filter to interpolate the original data onto a regular grid (section 2) and linear interpolate to fill gaps in the resulting interpolation. To reduce aliasing effects, sections are detrended and a Blackman window is applied. To allow averaging of spectral densities and to increase frequency resolution, we have padded the sections with zeros to a length of 750 km prior to applying the FFT algorithm.
We explored an alternative method for obtaining spectral densities, that of taking the Fourier transform from the biased autocorrelation function, but abandoned this approach in favor of the DFT. For complete time series, the autocorrelation method is equivalent to the DFT, but for data with missing values, it has the advantage that no interpolation is required across data gaps; hence, it appeared to be well suited to this project. However, tests with synthetic data (discussed in appendix B) indicate that the DFT gives more precise estimates than the autocorrelation function method for time series similar to those in our datasets.
2) Structure functions
We also estimate the second-order structure functions for the temperature sections. Structure functions provide a different perspective on the data than spectral density. In addition, it may be possible to estimate the slope of the spectral density using (1), when the structure of the data (i.e., the number and coherence of missing data points) makes it difficult (impossible) to employ standard techniques for estimating spectral densities directly (e.g., McCaffrey et al. 2015). An advantage of structure functions is that it does not require interpolation over missing data points; hence, even time series with relatively few data points improve the estimate of the average structure function. Ensemble-averaged structure functions are calculated by averaging over the contributions from all data-point pairs in the dataset for a given point separation (rather than by averaging over the structure functions computed for each section).
3) Slopes and confidence intervals
The slope of the spectral density or structure function over a given interval is estimated by fitting a first-order polynomial to the of the spectral density or structure function versus the of the wavenumber or separation of data points over that interval. The intervals used are noted in the text and in Tables 1 and 2. We estimate 95% confidence intervals by using a bootstrapping method (e.g., Efron and Tibshirani 1986) to generate 104 realizations of the average spectral densities.
In this section, we report the results of the analysis of the four datasets described in section 2. We start by comparing the three in situ datasets from the Oleander.
a. TSG, TEX, and ADCP
1) Average spectra
Figure 3 shows the DFT spectral densities of the ADCP, TSG, and TEX sections in the Sargasso Sea and Gulf Stream regions. A very strict quality standard of is used for all in situ data due to the large number of high-quality sections available. With that choice the number of sections used is 522 and 527 for ADCP, 205 and 178 for TSG, and 236 and 236 for TEX in the Sargasso Sea and the Gulf Stream, respectively. Spectral slopes are estimated for wavelengths between 10 and 100 km [0.1–0.01 cycles per kilometer (cpkm) ]for all three datasets, and for the interval 1–10 km (1–0.1 cpkm) for TSG and TEX, and are reported in Table 1. Generally, the Gulf Stream region is more energetic than the Sargasso Sea, and the slopes are steeper. Samelson and Paulson (1988) and Wang et al. (2010) suggest that spectral densities of temperature and potential energy are proportional. If that assumption is valid, then the slopes suggest that the Gulf Stream is somewhat closer to a quasigeostrophic regime, where a slope of would be observed, whereas other processes are of relatively larger importance in the Sargasso Sea. On the other hand, it has been shown that the slopes of density (potential energy) and temperature spectra can diverge due to density compensation (e.g., Rudnick and Ferrari 1999; Kolodziejczyk et al. 2015).
Although the general shape of the spectra is similar for the three datasets, their details differ. The ADCP sections are less energetic than the TEX and TSG sections, particularly at smaller scales ( km); thus, the ADCP spectra have steeper slopes. Qualitatively, this discrepancy can be explained by a response time of 5–10 min for the ADCP thermistor ( appendix A; Fig. 3), which is expected based on the location of the temperature sensor within the ADCP. The simple model [(A1)] used for simulating the ADCP thermistor response overestimates the damping of small scales relative to that of larger scales, suggesting that the thermistor response is more complex than assumed by (A1). Other issues may contribute to the discrepancy, albeit likely less strongly: The data products made available involve different filtering and averaging, and the ADCP temperatures are measured in the sea chest, whereas the thermosalinograph samples from the engine intake.
Spectral densities for TSG and TEX are consistent over a wide range of scales, with TEX being more energetic at scales of and smaller. Because TEX samples are taken directly at the water intake, modifications to the temperature prior to the measurements are expected to be minimal. This may not be true for TSG temperatures, which are on average biased 0.24 K warm compared with average TEX temperatures, suggesting that water warms within the system. This imposed heat flux not only increases temperature but likely also dampens short-term fluctuations, which reduces spectral energy at smaller scales. Despite these modifications, the TSG and TEX spectra compare relatively well, suggesting that the internal TSG sensor on the Oleander generates a robust measurement of the temperature structure, providing useful information when TEX records are unavailable.
2) Structure functions
The structure functions computed for the same section as the spectral densities in Fig. 3 are shown in Fig. 4. In general, a picture similar to that associated with the spectral densities emerges, with the Gulf Stream region being more energetic and showing steeper slopes than the Sargasso Sea. The structure functions also reproduce the differences among the three datasets seen in the comparison of their spectral densities; that is, TSG and TEX are more energetic at smaller scales.
Because (1) applies only when the spectral slope is constant and extends over an infinite frequency interval, it is not straightforward to compare slopes of the structure function at different separation intervals to spectral slopes at different wavelength. It has been suggested, however, that different high- and low-frequency regimes can be compared to small- and large-separation parts of the structure function if both parts are sufficiently well resolved (McCaffrey et al. 2015). Here, we explore that possibility by computing the slopes of the structure function λ over the intervals from 5 to 50 km and from 0.5 to 5 km (Table 2). Spectral and structure function slopes are compared graphically in Fig. 5. Generally, steeper structure functions correspond to steeper spectra, and changes in slopes from small to larger separations are consistent with those from small to large wavelengths, consistent with McCaffrey et al. (2015). It is also apparent, however, that λ is significantly flatter () than would be expected from (1), consistent with findings based on the simulated time series reported in appendix B. Furthermore, the steeper the slope the greater the discrepancy; that is, the discrepancy is generally more pronounced when the fits are over shorter spatial scales (1–10 km) and in the Gulf Stream region. This is consistent with the results of Huang et al. (2010), who show that the structure function is not a good method for estimating spectral slopes when there is a large difference in energy between small and large scales, that is, a steep spectral slope.
3) Seasonal cycle
An understanding of the seasonal variability in the spectral slopes is important because the VIIRS data distribution is not uniform over the year. Figure 6 shows the TSG spectral densities averaged by season over the dataset. Interestingly, the seasonal cycles in the Sargasso Sea and Gulf Stream regions have a different, almost opposite phase. In the Sargasso Sea, spectral energy reaches a maximum in summer. The relatively large seasonal differences, in particular at smaller scales, may explain part of the increase in the confidence interval in the mean spectrum at smaller scales. This is consistent with the findings of Luce and Rossby (2008), who used the 150-kHz Oleander data to estimate the seasonal variability of the eddy population in the Sargasso Sea and found that eddy numbers peaked in summer and fall.
On the other hand, Callies et al. (2015) found that potential energy spectra computed from data from the LatMix experiment in the Sargasso Sea are more energetic in winter than in summer (they also compute kinetic energy spectra from LatMix data and from the Oleander ADCP 75-kHz data at 50-m depth, which are also more energetic in winter.) Several issues may contribute to this discrepancy: First, the seasonal cycle of near-surface potential energy may be dominated by changes in mixed layer depth and stratification, rather than surface density variability. Furthermore, density variability may also be affected by density compensation (e.g., Rudnick and Ferrari 1999; Kolodziejczyk et al. 2015). Hence, it is not clear how well the surface temperature variability represents potential energy.
Second, the winter and summer LatMix experiments occurred in different areas; that is, the observed difference in seasonal cycles could be dominated by regional differences. Both LatMix areas are also closer to the Gulf Stream than our Sargasso Sea region, and in the Gulf Stream region we find that spectral energy is maximal in winter (right panel of Fig. 6; the blue curve lies just under the green curve there). Third, the surface quasigeostrophic mode decays exponentially from the surface, and its impact may be far smaller at 50 m (Callies and Ferrari 2013), where the kinetic energy spectra from Oleander data are computed, than at 5–6-m depth, where TEX temperatures are measured. In summary, our findings suggest that the seasonal cycle of submesoscale processes varies spatially, between different dynamical regimes.
4) Effects of gaps
To test the impact of gaps in real data, we have artificially introduced gaps into the TSG sections in the same manner as described in appendix B and recalculated the spectra and spectral slopes. Figure 7 shows 2D histograms indicating the change in slope resulting from gaps in the Sargasso Sea and Gulf Stream regions. As with the synthetic data ( appendix B), the spectral slopes are increasingly biased low as Q decreases and C increases, and this effect becomes more pronounced as the true spectral slope increases. These results motivated our choice to include only VIIRS sections with (dashed lines in Figs. 7a,b) in the following section.
1) Spectral densities
The spectral densities of VIIRS sections with and are compared to those of TEX in Fig. 8. With these criteria, the number of sections was limited to 65 (64) in the Sargasso Sea (Gulf Stream) region. The seasonal distribution of VIIRS sections in the Sargasso Sea is significantly skewed toward summer: The 65 sections available include only 2 in December–February (DJF), but 24 in March–May (MAM), 25 in June–August (JJA), and 14 in September–November (SON). Because Fig. 6a indicates a strong seasonal cycle in the Sargasso Sea, with spectral energy peaking in summer, particularly at smaller scales, this sampling bias is expected to affect the results. For this reason, we compare the VIIRS spectra to the TEX spectra in JJA (red curve in Fig. 8a; the annual mean TEX spectrum is also shown for reference). With the seasonal bias removed, we find good agreement on larger scales; however, VIIRS spectra are significantly more energetic on scales 20 km.
In the Gulf Stream region, there are also more sections in summer than in winter (18 vs 10). Since the distribution is more uniform and seasonal differences are weaker (Fig. 6b), we compare annual mean spectra in the Gulf Stream region (Fig. 8b). The spectra agree at scales 5 km, but the VIIRS sections are more energetic at smaller scales. Interestingly, although TEX and VIIRS spectra diverge at different spatial scales in the Sargasso Sea and the Gulf Stream, they diverge at the same energy levels (cf. Figs. 8a,b). Note also that because of the higher energy at smaller scales, VIIRS spectral slopes are flatter than those obtained from TEX (see Table 1).
Because VIIRS measures the ocean skin temperature (the upper 10 μm or so of the water column), whereas TEX measurements are taken at 5–6-m depth, the flatter and more energetic VIIRS spectra may be due to different physical signals [it is tempting to suggest that surface quasigeostrophic processes (; Blumen 1978) decay rapidly enough with depth so that the deeper TEX measurements reflect mainly quasigeostrophic processes (). The presence of a diurnal cycle in the upper 5–6 m, with homogenization during the night and restratification during the day, makes this argument problematic because it is difficult to reconcile the surface quasigeostrophic adjustment time scales with the diurnal cycle]. On the other hand, the VIIRS spectra could also be affected by the less stringent Q–C standards used, the higher noise level of the VIIRS sensors, or flaws in the retrieval process, which are discussed next.
To illustrate the expected effect of gaps on the VIIRS spectra, we have subsampled the TEX sections (only JJA in the Sargasso Sea) to the VIIRS resolution and introduced artificial gaps with Q and C sampled from the corresponding VIIRS sections with and . The spectral densities estimated from these sections are indicated by the green curves in Fig. 8. In the Sargasso Sea, the results suggest that the gaps explain only a relatively small fraction of the difference between VIIRS and TEX at scales smaller than 2 km (the green curve deviates from the red curve only at small scales in Fig. 8a). In the Gulf Stream region, on the other hand, the gaps explain a larger fraction of the discrepancy, as is expected because of the steeper slope of the spectral density. We repeated the calculation with stricter Q–C criteria to reduce the possible impact of missing data; however, the results are not conclusive because the significant reduction in the number of usable sections has the effect of increasing the confidence interval.
To illustrate the effect of noise on the spectral densities, we added white noise of amplitude 0.1 K7 in addition to the artificial gaps to the subsampled TEX sections (magenta curves in Fig. 8). Because the effect of white noise on the spectral density decreases as the ratio of signal to noise increases, the impact in the energetic Gulf Stream region is small compared with that in the Sargasso Sea, where it leads to a significant flattening of the spectra at smaller scales.
Another difference between the satellite and in situ data is that the satellite samples the entire area within a few minutes, whereas the Oleander requires about 1.5 days for one section. The impact of this difference is likely minimal owing to the rapid ship speed of 8 m s−1; the Oleander crosses 30-km mesoscale features in 1 h.
In summary, we note that TEX and VIIRS spectra diverge at different spatial scales but at the same energy levels in the Sargasso Sea and the Gulf Stream. This is a strong indication that sensor and retrieval noise is underlying these discrepancies. On the other hand, the satellite measures the ocean skin temperature, whereas the thermosalinograph samples at 5–6 m below the surface. Hence, it is possible that processes underlying the temperature fields differ at these depths.
2) Structure function
Structure functions for the VIIRS data are shown in Fig. 9. Since the structure function (7) is insensitive to gaps, all sections with are used, which increases the number of used sections to more than 150 in both regions. The structure functions show the same relations between VIIRS and TEX as the spectral densities. We further note that again, the slope of the structure functions reported in Table 2 is smaller than predicted by (1).
We also show the structure functions using only data from daytime [71 (74) sections in the Sargasso Sea (Gulf Stream)] and nighttime [83 (84) sections in the Sargasso Sea (Gulf Stream)]. We find that although the shape does not differ much between day and night, the temperature field is more energetic at daytime at all scales. This difference results from a difference in either the underlying daytime and nighttime SST fields or the processing. In cases with high diurnal warming, such as the eastern Mediterranean during summer months, we have noticed a significant increase in the mean SST gradient magnitude due to patchiness in the surface wind speed. This explanation is more likely to be valid in the Sargasso Sea, where the day–night difference is more pronounced in summer (not shown), as is the diurnal warming (Cornillon and Stramma 1985). It is not so likely in the Gulf Stream, where the day–night difference is independent of season. An alternate explanation is that the difference results from the different form of the VIIRS SST algorithm used for nighttime and daytime retrievals. Nighttime SST retrievals are based on three spectral channels—two in the 10–12-μm atmospheric window (M15 and M16, the VIIRS band or channel specifications) and the 3.7-μm band (M12)—while daytime retrievals rely only on M15 and M16; M12 was excluded to avoid sun glint in this spectral band. The difference in the set of spectral channels used in the retrieval process also results in differences in the portion of the algorithm that flags bad data—presumably pixels contaminated by clouds. It is not clear at this time which of these is the primary cause for the day–night difference.
In the present study, we calculate spectral densities and structure functions of surface temperature along sections between New Jersey and Bermuda. We use three in situ datasets: temperature measurements from an ADCP, a TSG, and its TEX. All instruments are mounted on the Oleander, a container ship making weekly round trips between New Jersey and Bermuda. In addition, we analyze satellite-derived SST from the VIIRS carried on the Suomi-NPP spacecraft. We compare two methodologies for computing spectral slopes from the datasets as well as with synthetic data, to estimate the effect of gaps, resulting from cloud cover in satellite-derived SST fields, on spectral slopes. These analyses, discussed in more detail in subsequent subsections, suggest that for SST in the Sargasso Sea and in the vicinity of the Gulf Stream,
the Oleander-mounted TEX and TSG provide consistent SST spectra down to scales of
the structure function is not a reliable estimator of the spectral slope in either region
spectra from the full-resolution VIIRS SST fields agree with in situ data on mesoscales but show elevated energy levels at smaller scales with possible physical and technical explanations
a. In situ data
Our analysis indicates that temperature spectra estimated from TEX and TSG compare well down to scales of . Spectra from TEX are more energetic at scales km, which we attribute to the warming of TSG temperatures while water flows through the intake pipes. Despite these drawbacks, our analysis suggests that Oleander TSG temperatures can be used as a reliable substitute for characterizing the temperature structure on larger scales (but not for absolute temperature) if TEX data are not available. Interior TSG data from other ships may also be useful; however, because of possible differences in the system (e.g., pipe length), the temperature modification may be different and damp the temperature structure at larger scales. Spectra estimated from the temperatures measured by the Oleander 75-kHz ADCP are significantly less energetic at submesoscales, in particular at scales km. This difference can be explained by the slow response time (5–10 min, appendix A) of the ADCP thermistor, which is due to the ADCP assembly.
b. VIIRS data
Spectra estimated from full-resolution VIIRS SST fields interpolated on a nominal Oleander track compare well with those from TEX for scales km in the Gulf Stream region and km in the Sargasso Sea region. At smaller scales, VIIRS spectra are more energetic and have flatter slopes. In both regions, the VIIRS spectra separate from the TEX spectra at a similar energy level, approximately 102 K2 km. This coincidence points toward pixel-to-pixel noise in the VIIRS SST fields as a possible explanation for the discrepancy, which is a typical problem for full-resolution SST fields derived from infrared satelliteborne sensors (Deschamps et al. 1981, 1984; Obenour 2013). On the other hand, simulated VIIRS spectra, obtained by adding data gaps and noise consistent with that in VIIRS to the TEX temperatures (based on pixel-to-pixel temperature differences from VIIRS fields, white noise with an amplitude of 0.1 K is assumed), cannot fully explain the difference between TEX and VIIRS spectra.
Since VIIRS measures the skin temperature of the ocean, whereas TEX samples at the engine intake at a depth of 5–6 m, it is also possible that different physical processes affect the temperature records of the two instruments. For example, the VIIRS data are more likely to be affected by diurnal warming. In this regard, we also found a significant difference between the daytime and nighttime temperature structures in VIIRS, in particular in the Sargasso Sea. Since it takes about 1.5 days for the Oleander to travel from New Jersey to Bermuda, the in situ sections include both, daytime and nighttime, the consequences of which are unclear. To better quantify the different causes for the difference between the in situ and VIIRS data requires further investigation.
c. Spectral density versus structure function
By comparing spectral slopes estimated by DFT and slopes of the structure function from our datasets, as well as from synthetic data with known properties ( appendix B), we have shown that (1) does not provide a good approximation of the relation between slopes of spectral density and those of the structure function for SST in the Gulf Stream and Sargasso Sea regions. Specifically, the structure functions were consistently flatter than expected from (1), and the error increases with the spectral slope. This is consistent with the findings of Huang et al. (2010): the structure function is not a good indicator of the spectral slope when significant amounts of energy are contained at the low-frequency end of the spectrum compared with the amount of energy at the high-frequency end.
An advantage of the structure function, on the other hand, is that it is not affected by data gaps; hence, for patchy datasets such as VIIRS SST fields for this region, the amount of usable data is significantly larger compared to DFT. For this reason, we found the structure function to be useful to contrast the structure of VIIRS temperature at daytime and nighttime and summer versus winter, which would have been prohibitive using spectra (e.g., in the Sargasso Sea, only two VIIRS sections satisfied the imposed Q–C standard in winter).
In conclusion, thermosalinograph and VIIRS data provide reliable spectra of the near-surface temperature of the ocean over a broad spectral range. However, several challenges remain, in particular at submesoscale and smaller scales for VIIRS. At such scales, the accuracy of the spectra is limited by sensor noise. This problem is more pronounced in the relatively quiet Sargasso Sea than the Gulf Stream due to the different signal-to-noise ratios. Further research is needed to better understand the spatial and temporal variations of submesoscale processes in the ocean, and, in particular, to resolve to what extent the difference at scales 20 km between VIIRS and thermosalinograph data can be attributed to noise versus differences in the temperature structure between the surface and 5-m depth. Despite these concerns, VIIRS is of significant value as a tool to study the structure of the upper ocean because of its global coverage—clouds permitting—from early 2012 to present. Finally, we note that spectral estimates may not be a good measure of short-scale features, such as fronts that cover only a small fraction of the ocean surface but that are readily resolved in satellite-derived SST fields, because of the sharp discontinuities with which they are often associated.
We acknowledge support by the Oleander Project and thank Tom Rossby, Charly Flagg, Alejandra Sanchez-Franks, Baylor Fox-Kemper, and two anonymous reviewers for their valuable feedback related to this work. We are grateful to NOAA AOML and the VIIRS project for making their data available. The work was also supported through NSF Grants OCE 0825845 and OCE 0851794. Salary support for P. C. was provided by the state of Rhode Island and Providence Plantations. B. B. was supported by Le Collège doctoral international de l’Université européenne de Bretagne.
Simulated ADCP Thermistor Response
To estimate the effect of the slow ADCP thermistor response, we simulate ADCP temperature assuming that TEX measures the true ocean temperature . We use a linear model,
where is the response time for the ADCP, with an initial condition of at .
To estimate a for the ADCP, we seek a best fit between simulated 3-min bin-averaged temperature and measured ADCP temperature by varying together with a constant shift in time and a constant temperature difference. Depending on the section and the cost function used, we obtain a between 5 and 10 min. This range exists because in comparison to the real ADCP, damping of high frequencies tends to be too strong relative to that of lower frequencies in the linear model (cf. Fig. 3), and the energy distribution between high and low frequencies varies between different sections and regions.
Tests with Synthetic Data
a. Synthetic data
We evaluate the performance of the methods described above by testing them with synthetic data with known spectral density. Specifically, we generate a synthetic time series by choosing a spectral slope γ and a length of the discrete time series N, with the number of contributing frequencies . To do this, we generate a random phase for each frequency and use the inverse discrete Fourier transform to generate a time series of length , of which we take the first N points. Varying the parameter allows us to explore the effect of energy below the resolved frequency interval on our spectral estimates. Gaps are generated by sampling from a lognormal distribution for a given ratio of good data points Q and coherence C, which determine the number of invalid data points, as well as the number of gaps in a section.
For these simulations we chose . Other parameter ranges are guided by typical values in our datasets: , typical of values in the open ocean, , , and . For each simulation, with a fixed γ, , Q, and C, we generate 100 time series, and the ensemble-averaged spectra and structure functions. Spectral and structure function slopes are determined from these averages.
For complete synthetic time series, the DFT generates the expected spectral densities with correct slopes. The method is relatively insensitive to the parameter ; that is, a variation of the number of frequencies used to generate the sample time series does not significantly alter the result, suggesting that application of the Blackman window minimizes aliasing.
Results differ for the structure function, as illustrated in Fig. B1a, which shows a 2D histogram that indicates the slopes of the structure function for separations from 5 to 256 grid elements for varying spectral slopes of the synthetic dataset and parameter . Generally, the slope of the structure function is less than predicted by relation (1); however, slopes increase and approach the prediction in the limit . This result is expected since in that limit (1) holds exactly. The difference between the estimated slope and the one predicted from (1) also increases with the spectral slope of the synthetic time series. This is because with larger spectral slopes, low frequencies become increasingly energetic relative to higher frequencies and hence have a greater impact on the variability. For the same reason, the error in estimating the slope also increases and structure function slopes are reduced as approaches N; that is, the energy of “unresolved” low frequencies is reduced (removed) from the signal.
When gaps are introduced in the synthetic data, the average structure function remains unchanged, as long as the number of data-point pairs available to estimate remains unchanged (not shown). The gaps, however, do induce a systematic error in the average spectral estimates obtained with the DFT. Figure B1b shows the slope of the spectral density as a function of a single gap of varying length in the middle of the synthetic time series and the “true” spectral slope γ. It illustrates that the estimate of the slope quickly degrades with increasing length of the gap, and more significantly with increasing γ. Next, we vary the coherence C of the missing data. Figure B2 shows 2D histograms of the difference between the slopes of the synthetic time series and those of the DFT estimate as a function of Q and C for . Again, estimates of average spectral slopes are precise for complete time series, but they deteriorate with decreasing Q. Figure B2 also indicates that estimates improve as C decreases for fixed Q, that is, as average gap lengths decrease.
Guided by these tests, we use a “Q–C standard,” which takes into consideration both, the fraction of good data points available (Q) and the structure of gaps (C), when selecting sections best suited for the spectral analysis of a given dataset. For each of the datasets considered, we choose a Q–C standard that is as strict as possible while retaining a reasonable number of sections to be used with DFT. Conversely, since the structure function is not sensitive to gaps, no Q–C standard needs to be applied.
A list of the acronyms used in this paper follows. (Additonal acronyms are available online at http://www.ametsoc.org/PubsAcronymList)
Acoustic Doppler current profiler
Atlantic Oceanographic and Meteorological Laboratory
Advanced Very High Resolution Radiometer
Discrete Fourier transform
Environmental data record
Fast Fourier transform
Joint Polar Satellite System
Scalable Lateral Mixing and Coherent Turbulence
Moderate Resolution Imaging Spectroradiometer
National Oceanic and Atmospheric Administration
National Polar-Orbiting Partnership
Ocean general circulation model
Ocean Surface Mixing, Ocean Submesoscale Interaction Study
Sea surface temperature………………………………………………………………
Suomi–National Polar-Orbiting Partnership
Remote temperature sensor
Visible Infrared Imager Radiometer Suite
The manufacturer estimates that the thermistor in the 75-kHz ADCP takes about 30–60 min to acclimate to a temperature change. The 150-kHz ADCP thermistor is believed to respond significantly faster, because it is glued directly to the conductive naval bronze housing.
At this point we are technically dealing with spatial series, but we continue to refer to them as time series, a nomenclature that is historically consistent.
We estimated an amplitude of ~0.1 K for the pixel-to-pixel noise in VIIRS SST retrievals from the distribution of pixel-to-pixel temperature differences in dynamically quiet regions.