Penetrative radiative flux in the Bay of Bengal

The Bay of Bengal (BoB), a semi-enclosed basin in the northern Indian Ocean, is a complex region with large freshwater inputs and strong vertical stratification that result in a shallow, spatially variable mixed layer. With the exception of shortwave insolation, the air-sea heat exchange occurs at the sea surface and is vertically redistributed by mixing and advection. Strongly stratified, shallow mixed layers inhibit vertical mixing, and the penetration of solar radiation through the base of the mixed layer can lead to redistribution of upper-ocean heat. This paper compiles observations of hyperspectral downwelling irradiance (Ed) from 67 profiles collected during six research cruises in the BoB that span a broad range of regions and seasons between 2009 and 2014. We report attenuation length scales computed using double and single exponential models and quantify the penetration of radiative flux below the mixed layer depth (Qpen. We then evaluate estimates of Qpen obtained from published chlorophyll-based models and compare them to our observations. We find that the largest penetrative heat flux (up to 40 of the incident Ed) occurs near 16Â°N where the mixed layers are shallow and the water is optically clear.

local warming, dynamical processes, and ocean-atmosphere interactions (Sweeney et al., 2005). Oceanic general circulation models typically parameterize upper-ocean radiant heating with two terms, one for the infrared wavelengths above 700 nm and the other for visible and ultraviolet wavelengths below 700 nm. Infrared radiation contributes 40-60% of the total incident irradiance and is absorbed primarily in the upper 2 m of the water column (Mobley, 1994). Ultraviolet and visible wavelengths can penetrate deeper and thus can play a more significant role in setting MLDs, overturning circulation, and low-latitude SSTs (Rochford et al., 2001;Murtugudde et al., 2002). Murtugudde et al. (2002) demonstrate that climate models can be improved by incorporating penetrative radiation and moreover that warming below the mixed layer by penetrative solar flux can weaken stratification and further enhance subsurface warming due to dynamic feedbacks resulting from weaker surface currents and reduced divergence. Shortwave radiant flux attenuates roughly exponentially with depth, and thus in regions like the Bay of Bengal (BoB), where oligotrophic waters and shallow mixed layers allow a significant fraction of energy to penetrate into the pycnocline, not accounting for penetrative solar flux may introduce substantial error when assessing mixed layer dynamics.
The BoB plays a significant role in the Indian monsoon and, in comparison with other ocean regions, it exhibits strong salinity stratification due to large freshwater input from major rivers such as Ganges-Brahmaputra, Krishna-Godavari, and Irawaddy (Sengupta et al., 2006). To represent shortwave heating in ocean models, it is essential to characterize the attenuation of penetrative solar radiation with depth. Variability in the attenuation of sunlight in the ocean water column primarily depends on the presence of optically active substances that absorb and scatter light, which are delivered through riverine discharge (e.g., sediment and dissolved organic matter) or produced locally within the euphotic zone via ecological processes (e.g., phytoplankton), in addition to the water itself. Here, we used five years of field observations of water column spectral irradiance, collected throughout the southern and northern BoB, to assess variability in the diffuse attenuation coefficient of downward irradiance at 490 nm (k d (490)) and of the photosynthetically available radiation (PAR) at wavelengths between 400 nm and 700 nm (k PAR ). Both coefficients provide insight into the attenuation of visible radiant flux in the surface ocean, necessary for refining our approaches to modeling mixed layer dynamics in the BoB. Such refinements are of considerable interest in ongoing efforts to obtain accurate SST predictions from models using satellite remote sensing.

INTRODUCTION
Ocean models of the tropical Indian Ocean often underestimate sea surface temperature (SST) tendency compared to observations, possibly due to the way in which atmospheric fluxes or ocean mixing processes are handled, or errors associated with the computation of penetrative heat flux (Chowdary et al., 2015;Ramu et al., 2016). The surface mixed layer SST tendency can be expressed as where, ∂T m /∂t is the rate of change of mixed layer temperature (T m ), ρ is seawater density (1,025 kg m -3 ), C p is the heat capacity of seawater (4,186 J kg -1 K -1 ), H m is the mixed layer depth (MLD), Q 0 is net surface radiant flux (W m -2 ), Q pen is shortwave radiant flux (W m -2 ) penetrating below the mixed layer, w e is entrainment rate (m s -1 ), and T b is the temperature at the bottom of the mixed layer. K z is the coefficient of vertical diffusion of heat (m 2 s -1 ). The divergence of the advective flux in the mixed layer is . (UT m ), where U is the horizontal velocity in the mixed layer, and z (~H m /2) is the depth between the mixed layer and the layer beneath. The vertical penetration of shortwave radiant flux through the upper ocean (Q pen ) has significant implications for " Here, we used five years of field observations of water column spectral irradiance, collected throughout the southern and northern BoB, to assess variability in the diffuse attenuation coefficient of downward irradiance at 490 nm and of the photosynthetically available radiation at wavelengths between 400 nm and 700 nm. " .
Hyperspectral irradiance profiles were collected using three different approaches on the six cruises ( Figure 2). On cruises SK-261, SN-30, SN-82, and SN-88, a profiling hyperspectral radiometer system (Satlantic HyperPro II) was deployed at each station to measure vertical profiles of upwelling radiance (L u (z, λ)) and downwelling irradiance (E d (z, λ)) calibrated at 136 wavelengths between 350 nm and 870 nm. The surface solar irradiance at the same wavelengths (E s (0 + , λ)) was measured with the system's above-water reference sensor mounted on deck. The HyperPro II was deployed away from the ship to avoid ship-induced perturbations and shading in irradiance (Mueller and Austin, 1995) and profiled with a vertical velocity between 0.4 m s -1 and 0.7 m s -1 , which provided a vertical resolution of approximately 10-12 samples per meter. Data were excluded when the profiler tilt was >5° to reduce any artifacts that would arise from the off-axis orientation during profiling. Raw instrument data were converted into physical units (e.g., W m -2 nm -1 for E d ) using a vendorprovided analysis package (ProSoft).
During the two remaining cruises, an E d (z, λ) sensor (TriOS Ramses ACC-VIS) was mounted either on the ship's profiling rosette (RR1317 cruise) or on a Wirewalker (RR1405 cruise). The Ramses ACC-VIS provides 190 effective channels between 305 nm and 950 nm, and an identical sensor was mounted to the upper deck of the ship to obtain FIGURE 2. Techniques for measuring the downwelling irradiance (E d ) varied between cruises. (a) On ORV Sagar Nidhi, measurements were made with a Satlantic HyperPro II system. This instrument falls slowly on a tether, as far from the ship and ship shadow as possible. (b) During the RR1 cruise, a TriOS Ramses ACC-VIS was mounted to the CTD rosette. Light interference by the ship required us to discard data from the upper 10 m. (c) During RR2, the TriOS Ramses was mounted on a Wirewalker (WW). The WW had a small surface footprint (~1 m) compared to the ship. Any single WW profile had variations likely due to wave focusing, clouds, and wobbling of the sensor during ascent, but when many profiles were averaged together, these errors were reduced and appeared reasonable within about 2 m of the surface.
Ocean Data View simultaneous measurements of E d (0 + , λ) at the same wavelengths. The Wirewalker (WW) is an autonomous vertical profiling system that uses wave power to ratchet downward, with no electronic components and only a few mechanical parts (Rainville and Pinkel, 2001;Pinkel et al., 2011;Lucas et al., 2016, in this issue). At the lower terminus of the profiling range, the cam that rectifies wave vertical motion is released, which permits the package to freely ascend, allowing high-resolution profiling of the upper water column. During RR1405, 1,025 E d (z, λ) profiles were collected using a Ramses radiometer mounted on a WW over a six-day drifting deployment (Figure 3). While single profiles sometimes exhibit significant variability, presumably due to steadiness of the sensor, wave focusing, and clouds, this effect was minimized by time averaging profiles over a period of time. Profiles collected between 10:00 and 14:00 local time were averaged to obtain daily E d (z, λ). Raw instrument data were converted into physical units using calibration data provided by the vendor and custom software written in MATLAB. For the subsequent analysis, all profile data were binned at 1 m vertical resolution and examined only within the euphotic depth, defined as the depth at which the PAR reaches 1% of its surface value. Data below this depth were excluded. In this study, the average depth of the euphotic zone at deeper-water stations was 65 (±8) m, whereas in the shelf region, sampled during SK-261, SN-88, and SN-82, it decreased to 54 (±15) m. The spectral downwelling irradiance, E d (z, λ), was extrapolated to the surface following a model defined in the SeaWiFS validation protocols (Mueller and Austin, 1995): where α represents the Fresnel reflection albedo for irradiance from the sun and sky (0.43), where 0 + and 0represent the depth slightly above and below the sea surface, respectively. These modeled E d (0 + , λ) were compared with those measured by the deck-mounted reference sensor that directly measured spectral incident irradiance E s (0 + , λ). Profiles that exhibited a difference between measured and modeled E s (0 + , λ) of more than 20% were discarded. This quality-control step resulted in a total of 67 profiles of E d (z, λ) available for subsequent analyses.

ATTENUATION COEFFICIENTS USING A DOUBLE EXPONENTIAL MODEL
The depth to which sunlight penetrates the oceanic water column is strongly dependent on wavelength. Red and near-infrared wavelengths are quickly absorbed by seawater and thus represent an important contribution to the surface heat budget, having a particularly strong influence on diurnal warming processes. Because these wavelengths are absorbed primarily within the upper few meters, it is necessary to have highquality irradiance measurements close to the surface to quantify their attenuation most robustly. The TriOS Ramses sensor measured irradiance over a large enough range (with factory calibrations up to 950 nm) to detect part of the nearinfrared waveband (~700-1,200 nm), whereas the sensors in the HyperPro profiler are calibrated from 400 nm to 700 nm only. Ship-shading effects during RR1317 precluded assessment of nearinfrared attenuation in these critical top few meters of the ocean, but use of the WW during RR1405 provided adequate profiling in the near surface to allow computation of composite daily profiles of spectral irradiance, resulting  best-fit estimates of ξ 2 and R (Table 1) were sorted and averaged depending on the season (winter or summer) and the latitude (above or below 15°N).

ATTENUATION COEFFICIENTS OF VISIBLE LIGHT USING A SINGLE EXPONENTIAL MODEL
To achieve as much consistency as possible between the different sensors and platforms, we examined light penetration in the BoB in terms of the PAR (the integrated downwelling irradiance from 400 nm to 700 nm) and the irradiance in the 490 nm band. Visible wavelengths penetrate more deeply into the water column than the near infrared, and the diffuse attenuation coefficients that are typically calculated to quantify their attenuation with depth (k PAR and k d (490)) are conceptually analogous to 1/ξ 2 from Equation 3 (Byun et al., 2014).
We computed k d (490) for all of the irradiance profiles as the linear slope of logarithmic E d (490) with depth, and computed PAR using The downwelling diffuse attenuation coefficient of PAR (k PAR )was then calculated using Equation 4 by replacing E d (z, λ) with PAR(z).
In our irradiance profiles, k d (490) varied between 0.03 m -1 and 0.15 m -1 and k PAR varied between 0.05 m -1 and 0.17 m -1 , with the highest values in the northern BoB shelf region (Figure 4). The northern BoB is greatly influenced by fresh water influx from perennial rivers such as the Ganges-Brahmaputra in the north, Mahanadi in the northwest, and Irawaddy in the east. Their influence extends up to the shelf region. This freshwater carries particles of organic matter and detritus that can absorb and/or scatter the visible light, resulting in a higher attenuation coefficient. The present study aims to use satellite data to estimate Q pen , which is a function of PAR and optical properties of the water column. PAR is in quantitatively adequate light profiles within a few meters of the surface. These curves were fit by a double exponential model (Paulson and Simpson, 1977) that represents the attenuation of long (nearinfrared) and short (visible) wavelength light separately using two attenuation length scales (ξ 1 and ξ 2 ): Because this model requires measurements in the infrared wavelengths, only spectral irradiance profiles from the WW deployments were used in this analysis. For each day of the WW deployment (n = 6), the optimal attenuation length scales (ξ 1 and ξ 2 ) and the relative magnitude factor (R) were calculated by minimizing the root mean squared error (RMSE) between the modeled profile (Equation 3) and the observations, using custom software written in MATLAB. Sea states and weather conditions were similar on all of these days, with overcast skies, no precipitation, and moderate wave conditions. The mean coefficients over this period, derived from this fitting analysis, are ξ 1 = 0.9 ± 0.4 m, ξ 2 = 20.8 ± 0.6 m, and R = 0.4.
The fitting procedure was then applied to all of the profiles from the Sagar Nidhi and Revelle cruises, keeping the parameters ξ 2 and R free, and fixing ξ 1 to the value determined from the WW data from RR1405 (0.9 ± 0.4 m). The resulting TABLE 1. Seasonal and regional values of the coefficients for the Paulson and Simpson (1977) double exponential model based on E d (z, λ) profiles measured from the Wirewalker. ξ 1 and ξ 2 are optimal attenuation length scales and R is the relative magnitude factor.  a standard ocean color satellite product (Frouin et al., 2012), and knowledge of k PAR is essential for estimating the propagation of PAR into the water column through the mixed layer. Therefore, we examined its relationship with k d (490), which is also a standard satellite ocean color product (Werdell and Bailey, 2005), and observed a high degree of correlation (R 2 = 0.89; k PAR = 0.0168 + 0.97 k d (490)) among the in situ data collected in the BoB for 67 profiles ( Figure 5). Presently, there are two operational standard algorithms for satellite estimation of k d (490) (Lee et al., 20015a). The first is a purely empirical algorithm that relates the blueto-green ratio of water leaving radiance to k d (490). In the second method, also empirical, the concentration of chlorophyll, a proxy of phytoplankton biomass, is first estimated using an empirical relationship based on the blue-to-green ratio of remote sensing reflectance, and the chlorophyll value is then used to estimate k d (490) based on another set of empirical relationships between k d (490) and chlorophyll. Semi-analytic algorithms have also been developed based on numerical solutions of radiative transfer, and provide an absolute percentage error of 14% in oceanic and coastal waters where k d (490) was in the range of 0.04-4.0 m -1 (Lee et al., 2005a). Therefore, we examined different algorithms for estimating Q pen , in the BoB to identify the most appropriate.

PENETRATIVE HEAT FLUX INTO THE BARRIER LAYER
The magnitude and distribution of heat in the water column affects fundamental upper-ocean processes such as air-sea fluxes and vertical and isopycnal mixing. In the BoB, solar energy that penetrates below the mixed layer into the barrier layer (the region immediately below the mixed layer) may become isolated from exchange with the surface due to stratification. Deeper mixing later in winter erodes the fresh mixed layer and mixes the sequestered heat upward. This process may influence cyclone tracks and intensity (Balaguru et al., 2012). In this study, we calculate the heat flux at the base of the mixed layer (Q pen ) from the mixed layer depth (MLD) measured at each E d (z, λ) profile. First, we computed a factor S ∫ 300 950 E d (0 + , λ)dλ PAR(0 + ) (6) S = from our estimates of E d (0+, λ) and PAR (0 + ) to quantify the relationship between incident PAR and the additional infrared radiation within the total incident radiation (Equation 6).
In order to arrive at spatially representative estimates of S, we used the surface reference data collected during RR1317, which covered a large transect across BoB. The results were very consistent from day-to-day, yielding S = 1.53 ± 0.04. Presuming that all of the near-infrared light is absorbed within the top few meters of the mixed layer, we predicted the fraction of penetrative heat flux relative to that incoming above the surface (Q pen /Q 0 ) as The MLD is defined here as the depth at which the density is 0.05 kg m -3 greater than at the surface. By this approach, estimates of the ratio of Q pen to Q o , ranged from 0.03 to 0.4, with the highest values tending to occur between 13°N and 17°N in the BoB, where mixed layers are shallow and fresh and the water is optically clear ( Figure 6).  We next computed the penetrative radiant flux (Q pen ) from published chlorophyll-based models and compared those results to our observations (Table 2).
We found that the Morel and Antoine (1994) and Ohlmann (2003) (MA94 and OH03) models were most correlated with the observations (R 2 = 0.73), and the intercept (0.17) was smallest in the case of OH03. Q pen computed using our parameterization showed mean ratio (0.86) and slope (0.83) closest to unity, with a comparable correlation coefficient (R 2 = 0.70), least RMSE (1.45%), APD (10.73%), RPD (7.75%), and UPD (6.89%) ( Table 3). The previous Q pen model parameterizations use chlorophyll as an input, considering it a proxy for phytoplankton attenuation of visible light. In the present case, we used an attenuation coefficient computed directly from observations. This fit reduces the error in the estimation of Q pen in the BoB, possibly because it encompasses the total contribution in visible light attenuation due to phytoplankton, detritus, and organic matter. TABLE 3. Statistical indicator obtained by comparing penetrative radiant flux at the mixed layer depth using present and published methods. The statistical indicator includes mean ratio; slope; intercept; regression coefficient (R 2 ); root-mean-square error (RMSE); and absolute (APD), relative (RPD), and unbiased (UPD) percentage difference between measured and modeled parameter.