The Nearshore Heat Budget: Effects of Stratification and Surfzone Dynamics

Temperature variability in the nearshore (from ≈ 6-m depth to the shoreline) is influenced by many processes including wave breaking and internal waves. A nearshore heat budget resolving these processes has not been considered. A 7-month experiment at the Scripps Institution of Oceanography Pier (shoreline to 6-m depth) measured temperature and surface and cross-shore heat fluxes to examine a nearshore heat budget with fine cross-shore spatial (≈ 20 m) and temporal (5 day to 4 hr) resolution. Winds, waves, air and water temperature, and in particular, pier end stratification varied considerably from late Fall to late Spring. The largest heat flux terms were shortwave solar radiation and baroclinic advective heat flux both varying on tidal time scales. The net heat flux is coherent and in phase with the nearshore heat content change at diurnal and semidiurnal frequencies. The binned mean heat budget has squared correlation R2 = 0.97 and best-fit slope of 0.76. Including an elevated breaking wave albedo parameterization reduced the residual heat flux and improved the best-fit slope. Baroclinic and barotropic advective heat fluxes have significant noise, and removing them from the heat budget improves the best-fit slope when stratification is weak. However, when daily mean stratification is large, baroclinic advective heat flux dominates variability and is required to capture large (≈ 3 ◦C h−1) internal wave events. At times, large heat budget residuals highlight neglected heat budget terms, pointing to surfzone alongshore advection of temperature anomalies.


Introduction
The nearshore region (defined here between the tidal shoreline and ≈ 6-m depth) includes the surfzone (where waves break) and adjacent shallow inner shelf. Nearshore temperature variability can reveal much about the dynamics and health of this important coastal habitat. For instance, temperature fluctuations can directly affect the nearshore ecosystem by altering growth rates, egg mass production rates (e.g., Broitman, Blanchette, & Gaines, 2005;Fischer & Thatje, 2008;Phillips, 2005.), coral health (Safaie et al., 2018;Schramek et al., 2018), and pathogen ecology and mortality (e.g., Goodwin et al., 2012;Surbeck, 2009). Temperature is a useful indicator of nutrient delivery to the nearshore (e.g., Omand et al., 2012), nearshore mixing and exchange due to (for example) rip currents (e.g., Hally-Rosendahl et al., 2014;Kumar & Feddersen, 2017) or internal waves, which can transport larvae and nutrients (Pineda, 1999) as far as the surfzone . However, the physical processes affecting nearshore temperature are not well-understood.
A range of physical processes can affect nearshore temperature. Advective processes such as internal waves can pulse cold (often nutrient rich) water into the nearshore region (e.g., Reid et al., 2019;Walter et al., 2014;Winant, 1974), and wind-driven upwelling can modify the thermocline depth and cross-shore location (e.g., Woodson et al., 2007;Walter et al., 2017). Air-sea heat flux (across the air-sea interface) includes shortwave and longwave radiant heat flux, as well as latent and sensible heat exchange which are often estimated with bulk parameterizations (e.g., Fairall, Bradley, Rogers, et al., 1996). Air-sea fluxes can establish diurnal cross-shore thermal exchange (e.g., Molina et al., 2014) and affect stratification and mixing (e.g., Price et al., 1986). Viscous dissipation of breaking wave energy is a unique surfzone heat source (Sinnett & Feddersen, 2014;Wei & Dalrymple, 2018), whereas breaking wave generated surfzone foam increases albedo (surface reflectivity) by as much as 8× (Sinnett & Feddersen, 2016) reducing shortwave radiation in the surfzone. The relative effects of wave heating and breaking wave albedo cooling depend on many factors including beach slope, wave height, latitude, season, and cloudiness . Surfzone breaking waves may also affect sensible heat flux (MacMahan et al., 2018) through Heat budgets relate heat content (temperature) variability with physical processes that transport heat into and out of the nearshore region. For example, on the Northern California coast, a heat budget in 60-to 130-m depths helped quantify how upwelling wind and wind relaxation events transport heat (Lentz, 1987). Similar observations demonstrated seasonal heat budget changes (Dever & Lentz, 1994). A heat budget for the New England Shelf (≈ 70-m depth) demonstrated how seasonal variation in air-sea heat flux and along-shelf advection affected temperature (Lentz et al., 2010). In shallower (≈ 12 m) depths, subtidal (30-hr period to 5 days) to seasonal heat budgets demonstrated summertime nearshore warming from surface gravity wave Stokes drift (Fewings & Lentz, 2011) and temperature response to stratification and alongshelf winds (Austin, 1999). Diurnal heat budgets have quantified wind forcing effects on a coastal embayment (Suanda et al., 2011;Walter et al., 2017). A heat budget analysis of a very shallow (≈ 1 m) reef environment demonstrated temperature sensitivity to depth and residence time (Davis et al., 2011) at hourly time scales.
However, the shallow nearshore often has large surface tides, internal waves (e.g., , surface gravity waves, and wave shoaling and breaking, all acting on subtidal, diurnal, semidiurnal, and higher frequency time scales. Thus a nearshore heat budget must include these time scales. The shallow nearshore is difficult to observe, and alongshore variation d∕dy is commonly ignored or estimated from sparse observations (e.g., Austin, 1999;Fewings & Lentz, 2011;Sinnett & Feddersen, 2014). A nearshore heat budget in 1-to 4-m depths closed on diurnal and longer time scales (Sinnett & Feddersen, 2014), although the temperature resolution was coarse, and advective processes were not resolved. Closing a nearshore heat budget is difficult, and consequently there are few shallow nearshore studies that involve both inner-shelf and surfzone processes. Many heat budgets ignore the very shallow inner-shelf and surfzone (e.g., Suanda et al., 2011), assuming these unobserved regions respond similarly to the deeper inner-shelf. However, cross-shore temperature gradients may interact with heat transport mechanisms (such as the nearshore internal wave field and rip currents) to modify exchange and transport dynamics Ignoring surfzone temperature gradients may affect resolution of transport and mixing dynamics extending onto the inner-shelf (Kumar & Feddersen, 2017;Grimes et al., 2019). This motivates development of a nearshore heat budget capable of resolving nearshore processes at time scales fast enough to capture internal wave variability.
Here, a nearshore (surfzone and inner-shelf) heat budget is tested using 7 months of observations (25 October 2014 to 1 June 2015) from a dense nearshore array, spanning depths from the shoreline to h = 6 m. Surfzone-specific wave heating and albedo modifications are included, and observations resolve internal wave driven changes to the nearshore heat content. The experiment and observational methods are outlined in section 2. The complete nearshore heat budget and terms are detailed in section 3. Observations of background conditions, individual heat flux terms, heat content variability, and unobserved heat budget residual are in section 4. Discussion of unique surfzone heat budget terms, internal wave runup, and the assumption of alongshore uniformity is in section 5. Section 6 is a summary.

Location
Nearshore heat budget observations were made during a roughly 7-month period (25 October 2014 to 1 June 2015) at the Scripps Institution of Oceanography (SIO) pier (32.867 • N, 117.257 • W). The 322-m long SIO pier extends west-northwest (288 • ) into water depth h ≈ 8 m (Figure 1a). A nearby canyon system bifurcates 1.7 km west of the pier-end, with Scripps Canyon extending to north-east and La Jolla canyon continuing to the south east.
Cross-shore bathymetry profiles at y = 0 m were taken along the SIO pier at 0.5-to 1-month intervals as conditions allowed. Near the SIO pier, the shoreline is approximately alongshore uniform with cross-shore slope ≈ 0.027 to depths near 20 m. The bathymetry changed throughout the observational period (see mean and standard deviation, Figure 1b) as waves and currents act to redistribute sediment (e.g., Ludka et al., 2015). Here, the nearshore domain extends cross-shore from the tidal shoreline at x sl along the pier into depth h ≈ 6 m at x os = −219 m, thus the full pier depth is not utilized. Mean tide level is the reference depth (z = 0) and the cross-shore coordinate extends positively toward the east along the pier with origin  (SIO pier)and surrounding nearshore waters with 10-m contoured bathymetry. The cross-shore (x) coordinate is aligned with the pier with origin at the mean tide level. (b) Detail of the cross-shore instrument deployment locations along the SIO pier (symbols) with reference to the mean sea level = 0 m (blue line), tidal standard deviation (blue dashed) and mean bathymetry (sand colored) with bathymetric variation (standard deviation, dotted). Grids (black lines) partition cells for heat budget estimation. Symbols indicate fixed locations of CTD (square), ADCP (triangle), TidbiT thermistors (blue circles), and SBE56 thermistors (red circles).
(x = 0) at the mean tide level shoreline. The alongshore coordinate origin (y = 0) is at the northern pier edge, extending positively toward the north (Figure 1b).

Data
Atmospheric and oceanographic observations were made throughout the nearshore domain. Air temperature, barometric pressure, wind speed and direction were reported at 6-min intervals from NOAA station 9410230 located on the pier near x os at elevation z ≈ 18 m. Humidity was observed at the pier-end in 5-min intervals (maintained by Earth Networks). A 4-way radiometer (Campbell Scientific NR01) mounted on a boom arm at x = −100 m extending 6.35 m to the south of the pier directly observed both downwelling and upwelling components of shortwave and longwave radiation at 1-min intervals (see Sinnett & Feddersen, 2016 for details).
Tidal elevation was recorded by NOAA station 9410230 as an average of 181 one-second samples reported at 6-min intervals with 2-cm accuracy and 0.1-cm resolution. A pressure sensor located at pier-end maintained by the Coastal Data Information Program (CDIP, station 73) reported hourly significant wave height H s and peak period T p . When these wave observations were not available (< 7% of the time), H s and T p were estimated from a spectral refraction model with very high skill initialized from offshore buoys (O'Reilly et al., 2016;O'Reilly & Guza, 1991. The pierend wave data and cross-shore bathymetry are used to model the cross-shore evolution of H s (x) and other wave properties such as wave dissipation and Stokes drift velocity U st . Salinity and water temperature were sampled every 6 min near the pier-end (x = −246 m and z = −5.8 m, square in Figure 1b) by a pier-mounted Seabird SBE 16plus SeaCAT operated by the Southern California Coastal Ocean Observing System. Salinity was linearly related to water temperature with weak salinity variation (< 1% variation 90% of the time). Nearshore temperature was observed by 33 thermistors (SBE56, and HoBo TidbiT) deployed on the SIO pier in the (x, z) plane (circles, Figure 1b). SBE56 thermistors sampled at 15-sec intervals and TidbiT thermistors sampled at 3-min intervals. All thermistors were calibrated before deployment in the SIO Hydraulics Laboratory temperature bath, obtaining 0.003 • C (SBE56) and 0.01 • C (TidbiT) accuracies. When thermistors failed (less than 1% of the time) the next nearest temperature observation at similar depth was used in its place. Current was observed by a bottom-mounted upward-looking 2 MHz Nortek Aquadopp ADCP located at x os but deployed ≈ 5-m north of the pier to avoid flow interference from pier pilings. The ADCP sampled with 1-min averages from z = −5.4 m (roughly 0.5 m above the bed) to below the surface wave trough in 0.3 m vertical bins. Velocity data were rotated into the x and y coordinate system.

Heat Content Changes
A nearshore heat budget is developed for a control volume bounded vertically by the bathymetry and wave-averaged free surface , and in the cross-shore by the fixed offshore location x os and the tidal shoreline x sl (Figure 1). The total nearshore heat content per meter of coastline is related to temperature such that, where the water density and the specific heat capacity of water c p are assumed constant, and the temporally and spatially varying temperature T is in Kelvin. The total nearshore heat content H tot (t) varies due to both changes in volume (mass) from tidal and sub-tidal sea level variation and changes in temperature.
This nearshore heat budget study is concerned with temperature-induced nearshore heat content variation (related to heat fluxes), rather than volume-induced heat content variation (primarily related to barotropic tides). Thus, temperature contributions are separated from tidal volume contributions by writing the heat content (1) as where T A is the average temperature in the domain and the cross-shore area of the domain is If domain temperature is uniform, the heat content time-derivative is where the 2nd term on the RHS is the volume-induced heat content change. As domain temperature variations are small (1% in Kelvin) relative to T A , the volume-induced heat content change is removed from further heat budget analysis by defining where dH∕dt represents the temperature (internal-energy) variability induced changes to the heat content of concern to this study.
Although volume changes are removed from the heat budget, they define a cross-shore tidal velocity U tide required to accurately account for advective heat flux contributions. Volume conservation inside the bounded alongshore uniform nearshore domain and a predominantly cross-shore tidal volume flux results in where the depth-uniform tidal velocity U tide is acting normal to the outer boundary at x os and the depth at x os is d x os (t) = h(x os , t) + (t). Thus, U tide depends on depth at x os and the tidal domain area so that

Heat Flux Balance
Nearshore heat content changes dH∕dt result from net heat flux F net through the sea surface and across the offshore and alongshore boundaries so that, where F air/sea is the heat flux exchange across the air-sea interface, F wave is wave energy flux, and F bt and F bc are the cross-shore barotropic and baroclinic advective heat flux components. For the average domain area (volume) in this study, a heat flux of 5.7 × 10 4 W m −1 would change the temperature by 0.1 • C in 1 hr. The heat flux exchange across the air-sea interface is, and consists of solar shortwave radiation Q sw , longwave radiation Q lw , latent Q lat , and sensible Q sen flux components. The onshore flux of mechanical wave energy F wave becomes a heat source due to the viscous dissipation of wave breaking (Sinnett & Feddersen, 2014;Wei & Dalrymple, 2018). Advective heat flux is decomposed into barotropic F bt and baroclinic F bc components (e.g., Austin, 1999).
Heat flux through the sea-bed is assumed negligible as geothermal and groundwater contributions are small, and the sand is well-fluidized. Here ∕ y of alongshore heat fluxes are assumed to be negligible. This assumption of alongshore uniformity is valid for the time-averaged alongshore momentum balance on beaches with alongshore uniform bathymetry and nearly incident wave field (Feddersen et al., 1998;Feddersen & Guza, 2003) and has been applied successfully to other heat budgets (e.g., Austin, 1999;Austin & Lentz, 1999;Sinnett & Feddersen, 2014). This assumption is discussed in section 5.3.

Heat Content Estimation
To estimate the nearshore heat content (1), the nearshore region is discretized into cells according to instrument location, with cell boundaries defined at the midpoint between thermistor locations (rectangles in Figure 1b). The heat content in each cell can then be estimated (e.g., Lentz, 1987) as, Here, the cell area A c is defined by the cell width x c = x 2 − x 1 and height z c = z 2 − z 1 , which is fixed for interior cells but is a function of (t) and h(t) for surface and bottom cells, respectively. The cell temperature T c is observed by each cell thermistor with = 1028 [kg m −3 ] and the specific heat capacity of water c p = 3993 [J kg −1 K −1 ] as derived for typical temperature and salinity from Fofonoff and Millard (1983). The total nearshore heat content per meter of coastline (1) is the summed heat content in all cells, The domain-averaged temperature T A required to remove tidal volume fluctuations as in (5) is an area weighted average of the thermistor array (see Figure 1),

Journal of Geophysical Research: Oceans
where Q sw d is the spatially uniform downwelling solar radiation component and is the cross-shore variable albedo (reflectivity) of the sea surface. The total cross-shore integrated heat flux contribution from shortwave solar radiation as in (9) is denoted here as, where ⟨·⟩ indicates cross-shore averaging and the nearshore width (ignoring wave setup effects) is The radiometer recorded 1 min Q sw d (t) averages which were then hourly averaged. Here, an albedo parameterization (Sinnett & Feddersen, 2016) based on solar zenith angle and breaking wave foam are used to estimate cross-shore (surfzone and the non wave-breaking inner-shelf) and time variable albedo. The parameterization utilizes a wave model and hourly averaged pier-end wave conditions (see section 2.2) to estimate foam induced albedo. The cross-shore averaging across the surfzone and inner-shelf (14) follows . As discussed in Sinnett and Feddersen (2016), with the weak observed wind speeds, whitecapping wave breaking albedo is negligible.

Longwave Heat Flux: Q lw
Net longwave radiation is where Q lw d is the time varying downwelling and Q lw u is the time and spatially variable upwelling components, respectively. Longwave radiation may be estimated for a blackbody as where is the sea surface emittance, T s (x, t) is the temperature of the ocean surface and is the Stefan-Boltzmann constant. The cross-shore integrated longwave radiation as in (9) is denoted here as, where ⟨·⟩ represents a cross-shore average. Here, ocean-entering (downwelling) longwave heat flux Q lw d was directly observed by the radiometer. Surface temperature T s (x, t) was estimated at each surface cell (see Figure 1b) and the emittance ≈ 0.98 following Josey et al. (1999).

Latent and Sensible Heat Flux: Q lat and Q sen
Cross-shore averaged latent Q lat and sensible Q sen heat flux are estimated from the COARE algorithm (Fairall, Bradley, Rogers, et al., 1996) with inputs including wind speed, air temperature, relative humidity, air pressure, cross-shore averaged surface water temperature, and salinity. Corrections for the cool-skin and warm-layer surface effects (Fairall, Bradley, Godfrey, et al., 1996) were made outside the turbulent and well-mixed surfzone. The cross-shore integrated latent heat flux contribution as in (9) is denoted Q lat L x . Spray droplets enhance latent heat exchange, however surfzone generated spray droplets are large and fall back to the ocean quickly and without exchanging a significant amount of latent heat (MacMahan et al., 2018). Here, potential surfzone breaking wave modifications to latent heat exchange are not considered.
However, large surfzone spray droplets may have an additional contribution to sensible heat exchange, dependent on the surfzone wave dissipation (MacMahan et al., 2018). Modifications to the surfzone sensible heat exchange Q sen were estimated following (MacMahan et al., 2018) and applied to the cross-shore averaged sensible heat flux, denoted Q sen L x . Surfzone sensible heat flux modifications to Q sen L x resulted in approximately a 3% increase above the COARE estimation.

Wave Energy Flux: F wave
The cross-shore wave energy flux due to viscous dissipation of surfzone breaking waves is a source of heat to the region (Sinnett & Feddersen, 2014;Wei & Dalrymple, 2018). For normally incident narrow-banded random waves, the cross-shore wave energy flux at x os is (Sinnett & Feddersen, 2014) where g is gravity, H s is the significant wave height at x os and c g is the group velocity. Wave reflection from shallow sloping beaches is small (typically < 3%; Elgar et al., 1994) and other pathways for dissipation of wave energy are either self-contained and frictionally balanced within the nearshore (e.g., wave-driven currents, bubble injection, sediment suspension) or negligibly small compared to breaking wave dissipation (e.g., sound and mechanical energy export from rip currents and undertow) . Thus, incident wave energy is assumed to completely dissipate as heat within the nearshore.

Advective Barotropic and Baroclinic Heat Flux 3.5.1. Velocity and Temperature at x os
At x os , the Lagrangian cross-shore velocity contributing to advective heat fluxes contains observed Eulerian U e and wave driven Stokes drift U st components (Fewings & Lentz, 2011), This cross-shore velocity at x os is then decomposed into barotropic and baroclinic components. As the tidally varying (volume) component is removed from the heat budget (5), so is the (barotropic) cross-shore tidal velocity (7) from the barotropic Lagrangian cross-shore velocity, making the barotropic velocity component Here, ⟨·⟩ represents a vertical average. Tidal U tide standard deviation was 0.2 cm s −1 (with zero mean) and vertically averaged Stokes drift mean and standard deviation U st = 0.6±0.5 cm s −1 . Lastly, the ≈ −2.0 cm s −1 long-term time average of vertically averaged U bt (likely from pier associated rip currents (e.g., Checkley & Lindegren, 2014) is removed to maintain nearshore continuity at long time scales. The baroclinic velocity component at x os is Since the ADCP vertical bin size was smaller than the cell heights (see Figure 1b), the baroclinic velocity in each cell U c bc is the vertical average of U bc (z, t) contained within that cell. Velocity observations at x os do not completely extend from −h to , as the ADCP requires a blanking distance above the transducer head and buffer from the turbulent and wavy sea surface and side lobe interference. Thermistors are at fixed vertical locations ( Figure 1b) and do not directly measure temperature at the bed or directly at the sea surface. Thus, T and U were extrapolated to include the entire vertical domain and allow for comparison to the heat content (which is also inclusive of the entire domain). Extrapolation maintained the local vertical gradient and was usually over short distances (over 1 m less than 10% of the time).

Cross-Shore Barotropic Heat Flux F bt
The nontidal barotropic cross-shore velocity at x os , U bt , is assumed to be balanced by alongshore gradients of the alongshore barotropic velocity V bt onshore of x os through volume conservation (e.g., Austin, 1999), resulting in where V bt is the spatially uniform alongshore component of velocity and A is the cross-sectional area. Thus, the estimated alongshore barotropic convergence required to maintain a constant volume is,

10.1029/2019JC015494
The net nontidal barotropic heat flux per meter of coastline includes flux contributions from the offshore boundary and the domain cross-sectional area (e.g., Austin, 1999), where ⟨·⟩ represents vertical averaging and T A is the average temperature in the domain. Combining (24) with (25) yields the barotropic heat flux F bt used in (8), where the boundary temperature difference is Alongshore velocity gradients required to maintain continuity are accounted for in (24), whereas alongshore T A gradients are assumed to be negligible.

Cross-Shore Baroclinic Flux F bc
The baroclinic cross-shore velocity U bc contributes a baroclinic heat flux at the offshore boundary, Here, the deviation from the vertical mean temperature at x os is T ′ x os = T x os − ⟨T x os ⟩. As velocity was observed at x os only, any alongshore baroclinic heat flux component is unobserved and reflected in the error. Vertical integration in (28) was estimated by summing the baroclinic contribution through each boundary cell at x os as, Here, T ′c x os is the deviation from the vertical mean temperature in each cell at x os and Δz c is the vertical height of each offshore cell.

Data Quality Control and Error
During the observational period, some data within the nearshore domain were either not collected, not valid, or outside the experimental assumptions required by our analysis. At these times, the data was excluded from further analysis. Thermistors deployed near the surface (z = −0.2 m) were frequently out of the water near low tide. When − H s ∕2 < −0.2 m wave troughs likely caused them to be exposed and their data was excluded. When thermistors at z = −0.2 m were submerged, surface cell temperature is the average of the two thermistors in the surface cell (see Figure 1b). Near surface velocity (ADCP) observations (above − 0.6 m) were excluded to prevent side-lobe interference. A small (≈ 5 degree) tilt error in the radiometer mounting axis was geometrically corrected using a time series of solar azimuth and zenith angles. When rain or heavy fog obscured the radiometer optics, shortwave Q sw and longwave Q lw heat fluxes were in error and were excluded (≈ 8 % of all observations). At these times and particularly during fall and winter, rain (frequently with elevated wind and cooler air temperatures) potentially contributed to relatively strong nearshore cooling, resulting in a mean dH∕dt bias.
Further data were excluded when extrapolated vertical temperature gradients (see section 3.5.1) were above 1.25 • C m −1 , increasing the potential for false signals (≈ 3 % of all observations). When alongshore velocity magnitude was above 0.12 m s −1 (two standard deviations from the mean), advective heat flux is more sensitive to small alongshore temperature gradients (which were not observed here). Data at these times were excluded (≈ 5 % of all observations).
To reduce the high frequency noise exacerbated by taking the derivatives in (5), hourly averaged heat budget terms are low-pass filtered with a 4-hr cutoff, zero-phase, 5th order Butterworth filter such that tidal variability is resolved. Filtered data is used where noted. Heat content error results from and bathymetry error (contributing error to A c ) and from thermistor error (see section 2.2). However, the A c error is much larger than T c error. Conservatively estimating a 7-cm bathymetric error (one quarter of the domain-averaged bathymetric standard deviation, Figure 1b) that varies on time-scales of weeks yields dH∕dt error of 10 3 W m −1 .

General Observations
Conditions observed at the experiment site varied on time scales from hours to seasons. Wave height H s (red, Figure 2a) varied between 0.2 and 2.2 m, with peak periods typically between 7 and 13 s (not shown). During this experiment, waves were uncharacteristically mild, with mean H s at the SIO CDIP station ≈ 10% lower than the 20-year average. The mixed barotropic tide included sea surface elevation fluctuations of ≈ ±1 m (black, Figure 2a). Tides combined with bathymetric variation caused the nearshore width L x to fluctuate between 132 and 256 m (Figure 2b). During large waves and low tides, the outer surfzone boundary, x sz , could be within 30 m of x os , though the domain always contained the entire surfzone.
Observed winds were typically calm, with hourly averaged wind speeds usually below 5 m s −1 (Figure 2c). Five storms (on 1 November, 12 December, 31 December -1 January, 27 February, and 25 April) elevated hourly averaged wind speeds above 8 m s −1 . Air temperature (black, Figure 2d) fluctuated between 29.2 • C in early November to 5.2 • C on 1 January. Hourly averaged surface water temperature at x os (red, Figure 2d) varied between 22.5 • C in late October to 14.0 • C briefly on 12 April. This variability was on seasonal, fortnightly, subtidal, diurnal, and semidiurnal time scales. In particular March-May had strong variability at fortnightly and shorter time scales. At x os , the surface temperature was warmer than the air temperature 78% of the time (compare black and red, Figure 2d). These observations were made during the warm "blob" (early 2014 to mid 2015) when coastal surface temperatures were uncharacteristically warm (Zaba & Rudnick, 2016). The small waves and weak winds at this site during the observational period are consistent with the larger scale "blob" conditions observed elsewhere in Southern California and the eastern pacific (Bond et al., 2015;Hartmann, 2015).
The hourly averaged vertical temperature gradient at x os between SBE56 thermistors located at z = −0.7 m and z = −4 m ( Figure 2e) had a strong difference in stratification, and thus the likelihood of internal wave generation, before and after 1 March. From 25 October to 1 March, dT∕dz at x os fluctuated between 0 • C m −1 and ≈ 0.4 • C m −1 . In contrast, after 1 March dT∕dz was often above 0.4 • C m −1 reaching as high as 1.5 • C m −1 (18 April, Figure 2e). The entire observational period is thus separated into a "Period I" before March 1st (containing generally low stratification) and a "Period II" after 1 March (often strong stratification).
The probability density function (pdf) of the vertical temperature gradient (dT∕dz) at both x os and farther onshore at x = −155 m quantifies the stratification increase between Period I and Period II, and the onshore reduction in stratification (Figure 3). Thermistors spanning the same vertical locations (z = −0.7 m and z = −2.5 m) at both x os and x = −155 m (red dots in Figure 1b) are used to estimate dT∕dz pdfs. In Period I, the dT∕dz at x os had a mean of 0.015 • C m −1 and a similar non-zero mode (Figure 3a) analogous to a Weibull distribution with shape parameter > 1. Farther (68 m) onshore at x = −151 m, the dT∕dz pdf (Figure 3c) had mean of 0.005 • C m −1 but the mode was zero (most often well-mixed), consistent with an exponential distribution. In Period II, at x os , the dT∕dz pdf shape was qualitatively consistent with Period I, but the mean and mode increased substantially to near 0.1 • C m −1 (Figure 3b). At x os , large stratification of dT∕dz = 0.4 • C m −1 was roughly 10× as likely to occur in Period II than Period I. During Period II, stratification also was reduced onshore. At x = −151 m (Figure 3d), the mean dT∕dz = 0.045 • C m −1 and the mode was Table 1 Mean near-zero (most often well-mixed) again consistent with an exponential distribution. The Period II dT∕dz pdf at x = −151 m was roughly 9× the mean in Period I, as the decay factor of the exponential distribution increased (Figure 3d). The dT∕dz reduction and change in pdf shape into shallow water reflects the stronger mixing and proximity to the vertically well-mixed surfzone.

Heat Flux Terms
The 4-hr filtered heat flux terms in (8) are examined over the entire observational period, and also separately for Periods I and II. Of the air-sea (F air/sea ) flux terms (9), shortwave solar radiation Q sw L x has the largest mean and largest variability (Table 1). Cross-shore averaged albedo ⟨ (x)⟩ varied between 0.1 and 0.43 and cross-shore integrated shortwave solar radiation ranged from 0 W m −1 at night to 1.81 × 10 5 W m −1 near noon on 18 May with mean and standard deviation Q sw L x = 2.97 ± 3.98 × 10 4 W m −1 (Table 1). Short-wave solar radiation largely varied diurnally, but also on fortnightly (due to clouds and albedo) and seasonal time scales. As available wintertime sunlight is lower than for spring, average Q sw L x was roughly 42% lower in Period I than in Period II (Table 1).
Other air-sea heat flux terms have a much smaller magnitude (Figure 4b, note change of scale). Cross-shore integrated longwave heat flux Q lw L x (magenta, Figure 4b) was nearly always negative, with mean and standard deviation −1.27 ± 0.56 × 10 4 W m −1 (Table 1), making it the largest net heat sink of all terms. Longwave heat flux had variability on fortnightly to diurnal time scales (Figure 4b) and longwave heat flux was similar in Periods I and II (Table 1). Cross-shore integrated latent Q lat L x and sensible Q sen L x heat flux (black and   (21) and the top to bottom difference in baroclinic velocity ΔU bc = U bc (x os , z top ) − U bc (x os , z bot ) (red), (b) boundary temperature deviation ΔT b as in (27) (black) and deviation from the vertical mean temperature as in (28)  green, respectively, Figure 4b) were also usually negative. Latent heat flux was the larger of the two terms, contributing −0.69 ± 0.67 × 10 4 W m −1 compared to sensible −0.15 ± 0.26 × 10 4 W m −1 (Table 1). Both latent and sensible heat fluxes had significant variability at fortnightly to diurnal time scales. The magnitude and variability of both Q lat L x and Q sen L x were larger in Period I when air temperature variability was slightly larger (see Figures 2c and 2d). Wave energy flux F wave (red, Figure 4b) is always positive and follows the wave height (see red, Figure 2a). During this unusually low-wave energy observational period, F wave contributed 0.25 ± 0.20 × 10 4 W m −1 to the cross-shore integrated surfzone heat budget (Table 1), the second smallest heat flux term. The F wave variability was largely at subtidal time scales. Usually at this location waves are larger in winter months, however, the F wave mean and variability was similar between Period I and II during this observational period.
Barotropic heat flux F bt (26) is composed of both velocity and temperature. Barotropic velocity U bt had a 0.02 m s −1 standard deviation (black, Figure 5a) over the observational period, an order of magnitude larger than U tide , with similar variability in Periods I and II. ΔT b (27) associated with F bt had a mean and standard deviation of −0.02 ± 0.10 • C (black, Figure 5b) with extremes of −0.7 • C and 0.5 • C. Characteristics of ΔT b are different in Period I and II. In Period II, the mean ΔT b value is 0.06 • C lower and the standard deviation is about 36% larger than in Period I. The resulting barotropic advective heat flux F bt (black, Figure 5c) had mean and standard deviation −0.34 ± 1.78 × 10 4 W m −1 . Although the mean contribution was relatively small compared with other terms, F bt had the third highest variability behind F bc and Q sw L x (Table 1).
Baroclinic heat flux F bc (29) is also composed of velocity and temperature. Baroclinic velocity U bc varies with depth and is illustrated by a vertical velocity difference between the top and bottom cells ΔU bc = U top bc −U bottom bc (red, Figure 5a). Over the observational period, ΔU bc = −0.007±0.037 m s −1 with extreme values above 0.1 m s −1 (Figure 5a). Unlike U bt , baroclinic U bc standard deviation is approximately 20% larger in Period II than in Period I, though baroclinic shear as high as ΔU bc ∕Δz = −0.072 s −1 can occur in both periods. At x os , the temperature deviation from the vertical mean depends on z, as illustrated by the difference between the top and bottom cells ΔT ′ = T ′ (x os , z top )−T ′ (x os , z bot ), (red, Figure 5b), similar to the vertical stratification (Figures 3a and 3b), with significantly more variation in Period II than in Period I (red, Figure 5b). Thus, both U bc and T ′ affect the resulting baroclinic advective heat flux F bc (red, Figure 5c), which had a mean and standard deviation 0.27 ± 3.67 × 10 4 W m −1 (Table 1).

Time scales of Shortwave and Advective Heat Fluxes
The total heat flux variability is dominated by solar shortwave (Q sw L x ) and advective (F bt and F bc ) heat fluxes ( Table 1). The Q sw L x spectra contains peaks at the diurnal and semidiurnal periods (black, Figures 6a and  6b), with ≈ 80% of the variance contained within the diurnal to semidiurnal frequency band. The Period II diurnal and semidiurnal solar spectral peaks are larger and slightly broader as springtime day length change  Figure 7. Spectra as in Figure 6 of total 4-hr filtered heat flux F net (red) and heat content change dH∕dt (black) as in (8)  was greater than the wintertime day length change. In Period I, F bt contains low frequency energy and a broad diurnal peak (black, Figure 6c) that was weaker than Q sw L x . During Period I, F bc (red line, Figure 6c) had similar variability to F bt , consistent with Table 1. During this low-stratification and weak internal wave period, F bc variability was mostly in diurnal to fortnightly frequency bands (Table 2). During Period II, the F bt standard deviation increased ≈ 50% (Table 1), which was also broadly distributed in frequency space (black, Figure 6d). The F bc advective heat flux increased most dramatically between Periods I and II, coincident with the increased stratification, with spectra elevated nearly an order of magnitude relative to F bt (Figure 6d). The Period II F bc spectra is largely white with diurnal peak. Although F bc variance has significant contribution at fortnightly and subtidal time scales, the largest increase is in the high-frequency band (Table 2), associated with nonlinear internal waves .

The Heat Budget
The net heat flux F net and nearshore heat content change dH∕dt are related through (8), (i.e., dH∕dt ≈ F net ) provided the terms are accurately estimated and neglected terms that involve the alongshore derivative ( ∕ y) are unimportant. The mean and standard deviation of both F net and dH∕dt are shown in Table 1, though the relationship between F net and dH∕dt is examined here in further detail. During Period I, spectra of nearshore heat content change dH∕dt and nearshore heat flux F net across time scales spanning 5 days to 4 hr have prominent peaks at diurnal and semidiurnal frequencies (Figure 7a). The dH∕dt and F net spectra are within the 95% confidence interval at all frequencies above 0.5 cpd (gray shading, Figure 7), indicating that F net and dH∕dt variability are similar across broad time scales. Diurnal and semidiurnal frequencies contain approximately 85% of the total dH∕dt variance in Period I. At these dominant frequencies, dH∕dt and F net are coherent above the 99% confidence interval, and the phases are zero with 90% confidence (e.g., Emery & Thomson, 2001).
In the more stratified Period II, dH∕dt and F net variability increased relative to Period I (compare Figure 7a and b), driven largely by F bc (Figure 6b). Most Period II variability was also at diurnal and semidiurnal frequencies (Figure 7b), albeit with broader peaks than in Period I, especially at the semidiurnal frequency. At these frequencies, dH∕dt and F net are again coherent and in phase. High frequency (above 4 cpd) dH∕dt variance was 16× stronger in Period II than in Period I and is discussed further in section 5.2.
Over the observational period, hourly F net and dH∕dt are reasonably correlated (R 2 = 0.48) with best-fit slope of m = 0.74 and intercept b = −0.44 × 10 4 W m −1 (gray dots, Figure 8a). For the entire observational period the binned heat budget (blue dots in Figure 8a) reduces noise and has binned mean R 2 = 0.97, best fit slope of 0.76 and intercept of −0.55 × 10 4 W m −1 (see Table 3). The best fit slope m < 1 indicates that estimated F net is biased high, or dH∕dt is biased low during extreme heat flux events. Bin standard deviation (vertical red lines, Figure 8a) indicates unresolved error variance due to estimation error or violated assumptions.
The different dynamics between Periods I and II are also reflected in the heat budget statistics. Both dH∕dt and F net have weaker variability in Period I compared to Period II as reflected by the spread of hourly values and bin centers (gray and blue dots, Figures 8b and 8c). During Period I, the binned mean best-fit slope is 0.94 (near one) and the binned standard deviation is small (vertical red lines, Figure 8b) indicating the Figure 8. Low-pass filtered heat budget terms: dH∕dt versus heat flux F net for (a) the entire observational period, (b) Period I, (c) Period II. Shown are hourly averages (gray dots), binned means ± standard deviation (blue dots and red lines), and linear binned-mean best fit (red dashed) with fit values in Table 3. Each bin includes 170 observations. The 1:1 line (black) is shown for reference.
heat budget closes when stratification is weak. In contrast, the hourly and binned-mean heat budget skill (Table 3) was lower in Period II than in Period I, reflecting the additional variability in both F net and dH∕dt at all frequencies. Still, the Period II binned mean heat budget is linear with R 2 = 0.99, with slope 0.70 suggesting an estimation bias (mean-error).

The Heat Budget Residual
The difference between the observed F net and dH∕dt is the heat budget residual (scatter of the gray dots around the 1:1 line in Figure 8) and is due to estimation error and violated assumptions such as neglecting ∕ y terms. Note, the residual is not related to times when both F net and dH∕dt were unobserved or excluded during rain events or when data were invalid (see section 3.6). Generally, temperature was accurately observed with excellent spatial and temporal resolution, although occasional inoperative thermistors resulted in temperature extrapolation to adjacent cells. Slowly-varing bathymetric error is estimated to lead to dH∕dt error of 10 3 W m −1 (see section 3.6), a few percent of the total observed variability. Thus, the heat budget residual is likely due to error in F net or from the assumption of alongshore uniformity. Definitively quantifying the source of the heat budget residual is impossible as alongshore temperature variations were not measured. Most (≈ 90%) F net variability is at diurnal and higher frequencies (see section 4.3). Longwave, latent, sensible and wave breaking heat flux may contain estimation error, though their overall magnitude Table 3 Heat

budget best-fit statistics: Squared correlation (R 2 ), slope (m) and intercept (b) for both the hourly data (left), and the binned means (right)
Hourly Bin Mean Obs Note. Rows separate statistics for the full observational period, Period I and Period II as shown in Figure 8.
is small and their contribution to variability at diurnal and higher frequencies is low compared to solar and advective flux.
However, advective (F bt and F bc ) heat flux spectra are relatively white and had large high frequency (< 1∕10 cph) variability, suggesting F bt and F bc may have relatively high noise levels. Estimating advective heat fluxes is challenging, particularly at tidal and higher frequencies and advection is frequently either filtered or averaged to isolate diurnal or tidal scales (e.g., Ulloa et al., 2018). Thus, the hypothesis that the advective heat fluxes are noisy is examined by reanalyzing the heat budget without F bt and F bc .
During Period I the heat budget improves if F bt and F bc are excluded from F net (compare Figure 9a with Figure 8b) with squared correlations (R 2 ) improving from 0.53 to 0.64 and best-fit slope closer to 1 (0.91). The binned mean heat budget similarly improves. Thus, during Period I, advective heat flux terms are noise dominated, likely due to ADCP noise from the typically very small observed currents and under-sampled near-bed and near-surface regions.
In Period II, the baroclinic signal is more energetic (Figure 6b) and the heat budget skill is higher if both F bc and F bt are included. Removing the advective terms reduces R 2 from 0.47 to 0.45, although the slope improves from 0.72 to 0.96. The binned mean standard deviation increases without advective terms (red lines, Figure 9b) as extreme values of dH∕dt are not well-described by F net . Thus, in addition to short-wave solar radiation, advective heat fluxes (especially F bc ) are important in Period II to capture the large dH∕dt events commonly associated with internal waves, as the larger baroclinic signal increases the signal to noise ratio. However, the reduction in best-fit slope by including F bt and F bc suggests that there may be an estimation bias. Further error may be due to violated assumptions as discussed in section 5.3.

Surfzone albedo and wave heating modifications
Solar zenith angle dependent parameterizations (e.g., Payne, 1972;Taylor et al., 1996) are commonly used and included in the Air-Sea MATLAB Toolbox (crusty.usgs.gov/sea-mat/), in the Finite Volume Community Ocean Model (FVCOM) (Chen et al., 2004), and in the Regional Ocean Modeling System (ROMS). Solar zenith angle derived albedo has typical values near = 0.06 (e.g., Payne, 1972;Taylor et al., 1996) which is frequently applied to the entire nearshore region when estimating solar radiation (e.g., Connolly & Lentz, 2014;Whitney et al., 2016). By also including breaking-wave albedo effects as in Sinnett and Feddersen (2016),  (see section 3.4.1), the estimated cross-shore averaged ⟨ ⟩ was elevated over standard (non-breaking) parameterizations by mean and standard deviation of 0.06 ± 0.03 and at times by as much as 0.19. Neglecting the breaking-wave surfzone albedo and instead applying a solar zenith angle derived albedo (Taylor et al., 1996) to the entire nearshore results in a best-fit slope 8% farther from unity and a 55% increase in heat budget residual. On average, applying a surfzone albedo correction reduces the residual magnitude by 0.38 × 10 4 W m −1 (not shown). Thus, accounting for the elevated surfzone albedo is critical to this nearshore heat budget where the surfzone constitutes a significant part of the area of interest. Further effects of cross-sectional width, waves, bathymetric slope, latitude, season, and clouds are discussed in . In regions with a highly reflective substrate or significant turbidity, further albedo corrections may also be warranted (Fogarty et al., 2018).
Including surfzone wave heating effects in this heat budget was less important. The F wave contribution is expected to be lower than breaking wave albedo effects in small wave conditions at mid-latitudes . During this experiment, waves were uncharacteristically small (section 4.1), reducing F wave . Small waves also reduce the width of the surfzone (where F wave acts), and thus the ratio of surfzone width to L x (on average ≈ 0.4) was generally smaller than usual. This resulted in smaller F wave contribution than in a previous experiment (Sinnett & Feddersen, 2014).

Internal Wave Related Baroclinic Advective Heat Flux
This location has frequent nonlinear internal waves (NLIW) that can propagate through the domain as bores to the surfzone . The NLIW bores bring near-bed cold water upslope in a manner consistent with an upslope gravity current . The significant F bc variability at semidiurnal and higher frequency time-scales suggest NLIW contribute to the heat budget and domain temperature variability. Although the 4-h filtered advective heat flux terms are inferred to have significant noise over the long durations of Periods I and II (Section 4.5), the 30-min estimated advective heat flux terms (section 2.2) can close a shorter duration heat budget during a Period II strong NLIW events. Here, an example internal wave event and the effect on the nearshore heat budget are explored using the high-temporal resolution observations.
An 18 hr period bracketing nighttime from 10-11 April 2015 highlights a typical Period II internal wave event (Figure 10). At this time, winds were weak (< 3 m s −1 ) and waves were small (H s ≈ 0.6 m) with strong spring tides (rather than the usual mixed tide). The 5 h from sunset to the large internal wave event at hour zero had weak cooling from longwave radiation (Q lw L x ≈ −1.5 × 10 4 W m −1 , magenta, Figure 10a) and both barotropic and baroclinic advective heat fluxes were relatively small (Figure 10b). Just before midnight, a large NLIW event propagated into the observational array. At x os , temperature at z = −1.4 m rapidly dropped ≈ 3 • C (red, Figure 10d) with isotherm displacements a significant portion of the entire water column, indicating a nonlinear internal wave. Concurrently, cross-shore current became strongly baroclinic with near bottom velocities ≈ 11 cm s −1 (onshore) and near surface velocity ≈ −10 cm s −1 (offshore) characteristic of a nonlinear onshore propagating bore. The event propagated onshore, cooling the nearshore over a ≈ 1.75 h "cooling phase" at a maximum of −1.4 × 10 6 W m −1 (black, Figure 10c). After an hour, the event Figure 10. Half-hour averaged time series of (a) shortwave Q sw L x (blue), longwave Q lw L x (magenta), latent Q lat L x (black), sensible Q sen L x (green) and wave energy F wave (red) heat flux, (b) barotropic F bt (black) and baroclinic F bc (red) heat flux, (c) the nearshore heat content time derivative dH∕dt (black) and net heat flux F net (red) and (d) associated temperature at z = −1.4 m in the surfzone (black) and at x os (red) for an 18-hr period surrounding an internal wave event beginning at 23:45 PDT on 10 April, 2015.
arrived at and rapidly cooled the surfzone at the same z = −1.4 m vertical level (black, Figure 10d). At 01:30 PDT, baroclinic velocity quickly reversed after the event reached its runup extent, consistent with earlier observations . At x os , the near bottom current switched from 5 to −10 cm s −1 (onshore to offshore) over 5 min, while near surface current switched from −5 to 7 cm s −1 (offshore to onshore). As the cool bottom water receded over the next ≈ 1.75 hr "warming phase", the nearshore warmed at a peak rate of 1.5 W m −1 (black, Figure 10c). The stratified inner-shelf warmed quicker than the well-mixed surfzone, so that the surfzone was at times as much as 1.5 • C cooler than the inner-shelf at the same near-surface z = −1.4 m level (Figure 10d). However, for both x os and the surfzone, temperature does not recover to its pre-event levels at 5-8 hr prior, consistent with the event described in . This indicates that the event is not reversible and can lead to longer time scale heat and temperature variability. Over the event's "cooling phase" and "warming phase," F net captured ≈ 80% of the observed reduction and gain in nearshore heat content, with F bc accounting for the great majority of this total and F bt accounting for the remainder (Figure 10b). After the event, dH∕dt variability was higher than before with notable dips and peaks around hours five to eight that were not captured by F net (Figure 10c).
With the stronger Period II stratification, such NLIW events (e.g., Figure 10) are common, causing F bc to be the most variable heat flux term while shortwave radiation Q sw L x is the second most variable term (Table 1). Both Q sw L x and F bc are most energetic at diurnal and shorter time scales (Figure 6), though over a typical daily cycle, the relative importance of Q sw L x variability to F bc variability is unclear. To investigate this relative importance, the standard deviations of Q sw L x , F bc and dH∕dt (denoted for example (Q sw L x )) are estimated Figure 11. Ratio of daily baroclinic heat flux standard deviation (F bc ) to shortwave solar radiation standard deviation (Q sw L x ) versus daily averaged vertical stratification at x os , dT x os ∕dz. The colorbar represents the daily heat content change standard deviation (dH∕dt). Standard deviation is estimated over a 24-hr period and each symbol (219 total) represents a day. Circles and triangles indicate Period I and Period II observations, respectively.
over 24-hr periods encompassing both semidiurnal and diurnal time scales. Daily mean stratification at x os , dT x os ∕dz, is also estimated as it is this stratification that NLIW propagate upon and that is modulated by shortwave solar radiation.
The (F bc )∕ (Q sw L x ) ratio depends on the stratification and affects the nearshore heat content variability ( Figure 11). Daily averaged stratification at x os , dT x os ∕dz, varied between well-mixed (10 −2 • C m −1 ) and strongly stratified (0.4 • C m −1 ). Note the 24-hr averaged dT x os ∕dz is smaller than the hourly averaged dT∕dz in Figures 2 and 3. When stratification is weak (as typical in Period I, circles Figure 11) internal wave events are rare and (Q sw L x ) is larger than (F bc ) by an order of magnitude. For these cases with dominant (Q sw L x ), dH∕dt variability [ (dH∕dt)] is relatively small (dark symbols, Figure 11). As stratification increases (typical in Period II, triangles Figure 11), F bc variability increases as well, so that when dT x os ∕dz ≈ 0.1 • C m −1 the relative variability of F bc and Q sw L x are equal. Above dT x os ∕dz ≈ 0.1 • C m −1 , (F bc ) usually dominates. When F bc variability dominates, (dH∕dt) is more likely to be large as well (lighter symbols, Figure 11). Note that when semidiurnal or higher frequency variability is excluded, the resulting daily (F bc ) is 4-10× smaller, particularly for the strongly Stratified Period II.
Over the observation period, daily (F bc ) was larger than daily (Q sw L x ) roughly a quarter of the time (Figure 11). However, in Period II (represented by triangles, Figure 11), (F bc )∕ (Q sw L x ) > 1 nearly half the time. Furthermore, F bc was larger than (Q sw L x ) in 23 of the 25 days containing the largest (dH∕dt). Thus, when daily averaged inner-shelf conditions are stratified above a critical level (here, roughly 0.1 • C m −1 ) offshore generated nonlinear internal wave events may have a dominant influence on the nearshore temperature variability on time scales from minutes to hours.

Potential Alongshore Non-Uniform Effects
In section 4.5, the heat budget residual was attributed to estimation noise, particularly in the cross-shore advective heat flux, and to neglected processes, particularly alongshore gradients in alongshore advective heat flux. Here, a large heat budget residual "event" is examined. The event occurred near 25 March 2015 03:00 PDT when rapid nearshore cooling was not captured by F net (compare black and red, Figure 12a). During the night and early morning of this event, F net is weak (−6 to +5 hr in Figure 12a). For the 9 hr prior to the peak cooling, winds were weak (≈ 4 m s −1 ) from 245 degrees and waves were small (H s ≈ 0.8 m) and essentially normally incident (2 degrees). No internal wave variability (as in Figure 10c) was observed at this location during this 22-hr span. Beginning near 02:00 PDT the nearshore began cooling rapidly, loosing heat at a peak rate of dH∕dt = −4 × 10 5 W m −1 though F net was very small (Figure 12a). This unexplained nearshore cooling rate was roughly 20% that of the internal wave event detailed in section 5.2. This event was also not reversible as both dH∕dt is almost entirely negative and post-event temperature was significantly cooler.
However, the nearshore did not cool uniformly in the cross-shore. The surfzone cooled at −1 • C h −1 near the event peak (gray box, Figure 12a and cross-section Figure 12b) and cooled 1 • C over the entire event.
In contrast, the inner-shelf temperature was largely unchanged or even warmed slightly during the event's peak (Figure 12b). Thus, the nearshore dH∕dt was primarily due to surfzone temperature variation. Together with essentially no F bt or F bc over this time period, this implies that the anomalous dH∕dt is caused by alongshore gradients in alongshore advective heat flux (e.g., y (vT)) within the surfzone, neglected in this heat budget (section 3.2). Non-uniform alongshore advective heat flux may contribute substantially to the heat budget residual and reduce heat budget skill. What drove this alongshore advection of a temperature anomaly is unknown. The coastline curves and there are two offshore canyons near the site (just outside of the view of Figure 1). Surfzone alongshore temperature variability of ±0.5 • C over 1-3 km has been observed (b) Cross-shore (x) and vertical (z) objectively mapped temperature difference ΔT associated with the dH∕dt over the 30-min period highlighted in gray in (a). Temperature observations were made at thermistor locations (black dots). The outer limit of surfzone wave breaking was at x ≈ −152 m (red dashed). Note in (a) the time series are not 4-hr low-pass filtered.
previously at similar sites. Inner-shelf alongshore temperature variability can be generated by tidal mixing fronts (Connolly & Kirincich, 2019) and obliquely propagating NLIW which can propagate all the way to the surfzone .

Summary
Observations from a dense sampling array along the Scripps Institution of Oceanography Pier (from 6-m depth to shore) were made for 7 months to examine a nearshore heat budget including the surfzone. Reduced solar shortwave radiation from breaking wave generated foam (and increased albedo) and wave heating from breaking wave viscous dissipation are surfzone specific heat flux modifications included in this heat budget.
Conditions were largely vertically well-mixed (dT∕dz ≈ 0.015 • C m −1 ) for the first ≈ 4 months of the experiment. During this time, shortwave solar radiation Q sw L x was the largest and most variable term, and diurnal and semidiurnal frequencies contained ≈ 85% of the dH∕dt variability. However, when stratification increased from March to June (dT∕dz ≈ 0.10 • C m −1 ) baroclinic heat flux F bc associated with internal waves became the most variable term, at times causing nearshore temperature to drop ≈ 3 • C in 1 hr before recovering. With higher stratification, dH∕dt variance at frequencies above 1∕6 cph was 16× stronger.
Over the observational period, hourly dH∕dt and F net observations are correlated (R 2 = 0.48) and are coherent and in phase at diurnal and semidiurnal frequencies. A binned mean heat budget reduces noise and has improved best fit statistics (R 2 = 0.97 and slope m = 0.76). As breaking wave induced foam elevates surfzone albedo, including a parameterized surfzone albedo correction improved the best fit heat budget slope and reduced the residual (unexplained) heat flux by 3.8 × 10 3 W m −1 . Noisy estimation of advective F bc and F bt contributed to the heat flux residual, and removing these terms from the heat budget improved the squared correlation and best fit slope when stratification was weak and F bc variability was low. However, when stratification was above ≈ 0.1 • C m −1 , F bc variability dominated and the squared correlation improved by including internal wave runup effects. The slope is farther from 1 when F bc is included, indicating a possible estimation bias. There are times when F net did not account for large dH∕dt variation, likely due to violation of the alongshore uniformity assumption.