Towards the Development of the National Ocean Service San Francisco Bay Operational Forecast System

The National Ocean Service (NOS), Center for Operational Products and Services installed a Physical Oceanographic Real Time System (PORTS) in San Francisco Bay during 1998 to provide water surface elevation, currents at PORTS prediction depth as well as near-surface temperature and salinity. To complement the PORTS, a new nowcast/forecast system (consistent with NOS procedures) has been constructed. This new nowcast/forecast system is based on the Finite Volume Coastal Ocean Model (FVCOM) using a computational domain, which extends from Rio Vista on the Sacramento River and Antioch on the San Joaquin River through Suisun and San Pablo Bays and Upper and Lower San Francisco Bay out onto the continental shelf. This paper presents the FVCOM setup, testing, and validation for tidal and hindcast scenarios. In addition, the San Francisco Bay Operational Forecast System (SFBOFS) setup within the NOS Coastal Ocean Model Framework (COMF) is discussed. The SFBOFS performance during a semi-operational nowcast/forecast test period is presented and the production webpage is also briefly introduced. FVCOM, the core of SFBOFS, has been found to run robustly during the test period. Amplitudes and epochs of the M2 S2, N2, K2, K1, O1, P1, and Q1 constituents from the model tide-only simulation scenario are very close to the observed values at all stations. OPEN ACCESS J. Mar. Sci. Eng. 2014, 2 248 NOS skill assessment and RMS errors of all variables indicate that most statistical parameters pass the assessment criteria, and the model predictions are in agreement with measurements for both hindcast and semi-operational nowcast/forecast scenarios.


Introduction
In 1998, the National Oceanic and Atmospheric Administration (NOAA)/National Ocean Service (NOS) installed a Physical Oceanographic Real Time System (PORTS) in San Francisco Bay to provide water surface elevation, current, near-surface and near-bottom temperature and salinity, and meteorological information to promote safe and efficient navigation in this area [1].As PORTS only supplies measured data at selected stations, NOS' National Operational Coastal Modeling Program (NOCMP) is developing an Operational Forecast System (OFS) to complement the service.
Many researchers have employed different versions of TRIM (Tidal, Residual, Intertidal, Mudflat Model) [2][3][4] for their purposes in the bay.For example, Cheng and Smith [5] employed the TRIM2D (two-dimensional TRIM) for the San Francisco Bay Marine Nowcast.The TRIM3D was recently applied by Gross et al. [6] to the entire San Francisco Bay.The UnTRIM, an unstructured version of TRIM3D, has also been applied to San Francisco Bay by MacWilliams and Cheng [7].
Two-and three-dimensional models have been extensively applied to investigate the hydrodynamic and morphologic processes in San Francisco Bay.report the existence of sand waves with heights in the order of 2 meters at the entrance of the Bay and consider coastal process evolution and the numerical prediction of severe storms on the coastline initially using the two-dimensional vertically integrated mode of the Delft3D-FLOW model [11].Uslu et al. [12] developed a very high resolution two-dimensional vertically integrated model for tsunami forecasts in this region.
Fringer et al. [13] developed the non-hydrostatic option SUNTANS (Stanford Unstructured Nonhydrostatic Terrain-following Adaptive Navier-Stokes) model, which has also been applied in San Francisco Bay by Chua and Fringer [14].
Some of these structured and unstructured models [15,16] have been employed to understand the role of stratification and baroclinic circulation on salt intrusion in the northeastern part of Figure 1.They focus on the dynamic interactions between fresh water from the Sacramento and San Joaquin Rivers and salt water from the open ocean.Their studies, as will be seen later in this paper, have great value in evaluating the advantage and disadvantage of using a flow or stage river boundary condition for these two rivers.
The primary objective of the NOCMP, however, is to develop and operate a national network of OFSs to support NOAA's mission goals and priorities.This ongoing San Francisco Bay Operational Forecast System (SFBOFS) will become a new member of the existing OFS family.Up to now, NOCMP has successfully developed CBOFS (Chesapeake Bay), DBOFS (Delaware Bay), TBOFS (Tampa Bay), NGOFS (Northern Gulf of Mexico), CREOFS (Columbia River Estuary) and other OFSs which along the Atlantic, Gulf of Mexico, Great Lake and Pacific Coasts.With the use of the COMF (Coastal Ocean Model Framework) on NOAA's High Performance Computer (HPC), each OFS automatically integrates NOAA's observing system's data streams and the forecast output from meteorological and basin scale ocean models to generate necessary model input forcings, and then perform hydrodynamic model predictions with such forcings.Also with COMF, these OFSs perform nowcast and short-term forecast predictions (48 hours in most case) of pertinent parameters which include water levels, currents, salinity, and temperature and disseminate them to users.A state-of-the-art numerical hydrodynamic model driven by real-time data and meteorological, oceanographic, and river flow (or stage) forecasts forms the core of the end-to-end system.For detailed information on the COMF refer to Zhang et al. [17].NOS CO-OPS is evolving to support two hydrodynamic models: ROMS for structured grid applications and FVCOM for unstructured grid applications.
As San Francisco Bay (Figure 1) has complex topography and shallow water features (the average water depth is less than 5 meters in the Bay), the well tested unstructured Finite Volume Coastal Ocean Model (FVCOM) [18][19][20] is employed as the core of the SFBOFS.Another reason to employ FVCOM lies in the fact that this model has already been ingested into the COMF-HPC as one of the major core models.Measured river flow data from USGS are used as river forcings for the five small rivers in blue.Measured river flow or river stage data at Rio Vista and Antioch may be used as river forcings, respectively, for Sacramento River and San Joaquin River.Technically, if river stage data are employed, the grid points across these two major rivers are treated the same way as those on the Pacific Open Boundary.This paper shows the major steps in how SFBOFS has been developed, assessed and put into quasi-operational status.First, FVCOM is briefly reviewed in Section 2, followed by an overview of the model's setup in Section 3 for the tide and hindcast cases.Section 4 presents the model's astronomical tide-only scenario simulation evaluation, while the hindcast skill assessments are described in Section 5.The COMF setup and assessment of the quasi-operational nowcast/forecast test are discussed in Section 6. Conclusions and discussion are given in Section 7.

The Model Overview, its Grid, and Subsequent Revisions
The physics of the FVCOM model and many aspects of the computational scheme are equivalent to the widely used Princeton Ocean Model.The FVCOM model solves the three-dimensional, vertically hydrostatic, free surface, turbulent averaged equations of motions for a variable density fluid.The model uses a triangular unstructured horizontal grid with a generalized sigma vertical coordinate.Dynamically coupled transport equations for turbulent kinetic energy, turbulent length scale, salinity and temperature are also solved.The two turbulence parameter transport equations implement the Mellor-Yamada level 2.5 turbulence closure scheme.The numerical scheme employed in FVCOM to solve the equations of motion is summarized in [18][19][20].
The FVCOM application to San Francisco Bay uses external forcing of water level, ocean density, wind, sea level atmospheric pressure, air temperature, relative humidity, downward long wave radiation, short wave radiation, and fresh water discharges entering the model domain.The model calculates water levels, three dimensional velocity, salinity, and temperature.
The horizontal grid structure and open boundaries are shown in Figure 2.This grid was developed using the Surface Water Modeling System (SMS) Version 10.1 as described by Brigham Young University Surface Modeling Laboratory and was based on the VDATUM grid developed for the coastal waters of North/Central California, Oregon and Western Washington [21].The open boundary of the San Francisco Bay grid was developed from this grid in the near shelf region external to the Bay.It was necessary to modify the VDATUM grid such that the outer boundary of the San Francisco Bay grid follows an approximate circular arc with one of the element sides near orthogonal to the boundary arc.The grid contains 102264 elements and 54120 nodes with a minimum depth of 0.2 m and maximum depth of 106.8 m [22].A uniform 20-layer sigma level vertical discretization was considered.
The following element quality checks were used: (1) minimum and maximum interior angles of 10 and 130 degrees, respectively, (2) maximum slope of 0.1, (3) maximum adjacent element area change ratio of 0.5, and (4) maximum number of elements connected to a node of 8.Note the slope corresponds to the maximum allowed gradient of the edge length inside the domain.The slope determines how fast the mesh size will increase toward the middle of the region.A small slope order of 0.1 means small meshes.The paving method was used, which uses an advancing front technique to fill the polygon with elements.Based on the vertex distribution on the boundaries, equilateral triangles were created on the interior to define a smaller interior polygon.Overlapping regions were removed and the process is repeated until the region is filled.Interior nodal locations are relaxed to create better quality elements.Several triangles were adjusted such that the minimum interior angle was at least 30 degrees to improve FVCOM stability.In addition, along the open boundaries, the element topology was adjusted such that each boundary element contained only one boundary side.The triangle lengths are sufficiently small, that a reasonable M 2 wavelength to grid size is obtained as shown in Figure 3, where element lengths decrease from 400 to 1700 m along the open ocean boundary to a near uniform resolution of order 150 m throughout the interior bays and into the lower delta.
Several modifications were made in the development of SFBOFS to Version FVCOM 3.1.6.It should be noted that if the HEATING_CALCULATED_ON options is selected then the AIR_PRESSURE_ON option must be selected.While the sea level atmospheric pressure field is needed for the heating calculations, its gradient does not need to be applied in the momentum equations.In fact, for tidal simulations this is not correct.For tidal simulations with the heat flux calculations selected, it is necessary to provide a constant sea level atmospheric pressure field (1013 mb).Also, if one selects AIRPRESSURE_ON = F in namelist, the flag FLAG_28 = −DAIR_PRESSURE in file make.incshould be noted.The bottom roughness fix reported by Warner [23] for wetting/drying was added in file brough.F.In model testing, with the min_depth as 0.05 m, the model ran successfully and works for the wetting/drying case in San Francisco Bay.A Newtonian damping sponge layer was implemented by Lettmann [24], which provides a more robust implementation of the clamped water level open ocean boundary condition.
In the shallow mud flat regions of the Bay, there was also an issue with overheating.As a result, subroutine vdif_ts.F was modified to limit the short wave radiation and total heat flux as a function of depth.For depths less than 10 m, the fluxes were set to zero.In this manner, the heat transfer is due to only advection and diffusion.There, the zeta1_eff and zeta2_eff parameters which control the attenuation of the short wave radiation are set never to be less than 30% of the water depth and therefore always allow attenuation.In total, the following routines are involved in the above modifications: 1. fvcom.F, mod_ncdio.F, mod_timeseries.F-air_pressure option or heating_calculated_on option.2. brough.F-bottom roughness with the Warner [23] wet/dry treatment.3. advave_edge_gcn.F, advave_edge_gcy.F, extuv_edge.F, mod_semi_implicit.F and vdif_uv.F-Lettmann [24] sponge boundary.4. vdif_ts.F and vdif_ts_gom.F-revised heat flux in shallow water.
The interaction between the hydrodynamic and the sediment-water interface, particularly in the shallow water mudflat areas, which occupy some 16% of the Bay surface area, is an area where further research is needed.Fang and Stefan [25] considered the dynamics of heat exchange between the sediment and the bottom boundary layer for several hypothetical lakes.They found that the direction of the heat transfer reverses frequently on daily timescales as well as following an overall seasonal cycle based on weather conditions at Minneapolis-St.Paul, MN.Smith [26] performed a series of heat budget studies in Indian River Lagoon, FL, to estimate the water-sediment heat exchanges using assumed values for conductivity and density.The study sought to characterize subseasonal heat fluxes and temperature changes in the sediment and overlying estuarine waters.
The bottom stress formulation in shallow water for wetting and drying has received continuing interest.Research by Xue and Due [27], Uchiyama [28], Oey [29,30], and Oey et al. [31] has indicated that the bottom drag coefficient must be adjusted if the water depth approaches the bottom roughness height.How to perform this adjustment is an area for further consideration.In the present version of FVCOM, the effective water depth used in the bottom friction formulation is limited to 3 m; e.g., when the actual water depth is less than 3 m, the depth used in the bottom friction formulation is set to 3 m.

Model Setup
Basically, the model needs reasonable specifications of the following four items to obtain skillful predictability.They are (1) River boundary forcing conditions, (2) open ocean boundary conditions, (3) initial conditions, and (4) surface forcings.Each of these model elements is discussed below.

River Boundary Forcing Condition Specification
There are seven rivers considered in the model.Traditional river discharge condition is used for the five small rivers that are not in the delta (the rivers with names in blue in Figures 2 and 3) area.These five rivers are the Petaluma River (2 m 3 /s), Alameda Creek (3 m 3 /s), Napa River (100 m 3 /s), Coyote Creek (2 m 3 /s) and Guadalupe River (3 m 3 /s) with approximate mean annual flows in parentheses.The two rivers in the delta are Sacramento River at Rio Vista (1000 m 3 /s) and San Joaquin River at Antioch (100 m 3 /s) with approximate mean annual flows in parentheses.Two different upstream boundary condition types were considered for the Sacramento and San Joaquin Rivers forcing specification.
In type one, the average daily flow of Sacramento River was used to specify the flow at Rio Vista (RIO), while the San Joaquin River flow was estimated as the total delta outflow (OUT) minus the Rio Vista flow (RIO).The measured data are from the California Department of Natural Resources' DAYFLOW [32].The average daily flows (note: negative flow indicates flow into the Delta from the Bay) are used only for the hindcast scenario testing.Minimum inflow and zero salinity were set up for the two rivers in low flow period when DAYFLOW's estimates may be suspect as noted by Oltmann [33].In type two, the water level surface elevations (stage) were specified at Rio Vista and Antioch for Sacramento River and San Joaquin River, respectively, similar to the previous work in this region by MacWilliams et al. [16].
Both flow and stage river boundary conditions were used in the hindcast for comparison purposes.Please see Schmalz [22] for details.However, in the nowcast/forecast system we only use a river stage forcing for these two rivers.This is because NOCMP's top priority is to support PORTS for navigation safety, and water level prediction is paramount.Previous studies and personal communication with Michael MacWilliams [34] have found that flow boundary condition may be more suitable if the focus is on salinity prediction.However, our major concern is surface water level as in PORTS, and the stage boundary condition is more suitable.

Open Ocean Boundary Condition Specification
The open ocean boundary of the grid (see Figures 2 and 3) is forced with a superposition of the subtidal water levels and predicted tides.The harmonic constants of M 2 , S 2 , N 2 , K 2 , K 1 , O 1 , P 1 , and Q 1 that are used to predict tide are derived from the Oregon State University Tidal Inversion Software (OTIS) for the West Coast (WC2010 1/30°) [35].
In the hindcast scenario, the subtidal water level signal at Point Reyes (see Figure 1 for its location) is used to prescribe the subtidal water level along the outer boundary.A revised sponge layer treatment at the open ocean boundary was considered.The salinity and temperature at the open boundaries were determined, with nudging, from NOAA's World Ocean Atlas 2001 [36].No velocities are prescribed along the open ocean boundary.
In the nowcast and forecast scenarios, subtidal water level open boundary conditions are generated from the NCEP's (National Centers for Environmental Prediction) G-RTOFS (Global Real-Time Ocean Forecast System) gridded operational products.The temperature, salinity and baroclinic current open boundary conditions are also generated from G-RTOFS.The most recently available products for the given time period are searched and used for adjustments of the open boundary conditions using the COMF-HPC.Several horizontal interpolation methods are implemented, and a linear method is used for vertical interpolation from G-RTOFS vertical coordinates to model vertical coordinates.Measured real time sea surface elevation data at Point Reyes are used for the subtidal water level adjustment along the open boundary.Similarly, measured temperature data at San Francisco (see Figure 1 for its location) are used for boundary temperature adjustment.The adjustment is the difference between the observation and the G-RTOFS prediction at the start of the nowcast/forecast cycle and is discussed in greater detail in Section 6.1.

Initial Condition Specification
For the 19-month hindcast initial condition specification, the salinity and temperature fields were developed for 1 April 1979 using the joint NOS and USGS historical circulation survey conductivity-temperature-depth (CTD) datasets and the model was started from rest.The quasi-operational nowcast/forecast system started in the middle of March 2013 when a climatological temperature and salinity file (with adjustment from observation) was used as the very first initial condition.For each nowcast/forecast cycle, the COMF-HPC will automatically find the most recent restart file as this cycle's initial condition (SFBOFS has four cycles a day).The details of the HPC-COMF can be found in Zhang et al. [17].

Surface Forcing Specification
For hindcast scenario, the North American Regional Reanalysis (NARR) 2007 datasets with 32 km spatial and 3 h temporal resolution were interpolated to the model grid to provide 10 m winds, sea level atmospheric pressure, and 2 m fluxes of downward shortwave radiation and net total heat flux.For the nowcast and forecast, the COMF-HPC will automatically find the most recent NAM4 (North American Mesoscale Model 4 km resolution) results in NCEP's data tank to get the necessary input surface forcings.

Tidal Simulation
The tide scenario simulation is the standard first step for all OFS' development.This is due in good measure to the fact that water level is the first priority for safe navigation, and tide and tidal current are the dominant dynamic processes in most coastal waters.For the tide scenario, the model setup for the four forcing specifications as mentioned in the previous section is similar to that for hindcast scenario.The slight differences can be found below.

Short Term Experiment: 1-15 April 1979
A three-dimensional simulation approach including baroclinics was used to capture the influence of internal waves on the tidal dynamics following [37,38].The slight model setup difference from hindcast scenario is: winds were set to zero and the sea level atmospheric pressure set to 1013 mb.River flow conditions are used for all rivers.The April 1979 NOS and USGS historical circulation survey data were used to compare the model results with the observation.
To develop initial salinity and temperature conditions on 1 April 1979 (and on 1 September 1980 for the later extended experiment case), the available CTD and CT time series data were placed on a coarse unstructured grid of order 50 elements.An interpolation program was developed in which each FVCOM grid node was assigned a given element and the salinity/temperature value interpolated from the node values at the appropriate depths.This program allows the initial density condition to be developed for the tidal and hindcast simulations.
To calibrate the bottom roughness, the approach of Cheng et al. [15] was used, in which the bottom roughness is made a function of water depth as in Table 1.To reduce the amplitude of the simulated water level response at Port Chicago, the bottom friction was further increased above Carquinez Strait as noted in Table 1.The water level response with respect to MLLW at Port Chicago for Experiments 1 and 2 is similar (See Figure 4).Results for Experiments 5 and 7 show very minor improvement in the agreement with water level observations at Port Chicago in the order of a 2 cm reduction in RMSE.Experiments 3, 4, and 6 were unstable, due to large horizontal gradients in bottom roughness during the wetting/drying cycle.Three additional Experiments 8-10 were conducted in which the river stage at Rio Vista and at Antioch was reconstructed from NOS harmonic constituents.Experiment 8 used the Experiment 7 bottom roughness specification.Experiment 9 included a 20 cm offset for the San Joaquin River and a 22 cm offset for the Sacramento River.In Experiment 10, the Experiment 9 offsets were retained and the Set 1 Bottom Roughness z 0 values were used.Note in these stage experiments the Oregon State University Tidal Data Inversion, OTIS Regional Tide Solutions [35] harmonic analysis results were reduced by 5% for the four ocean open boundary stations.Note the Sa and Ssa harmonic constituents derived from San Francisco water level analysis were used at these stations.All other open boundary node water levels were derived via linear interpolation of values from two of the stations surrounding the node.In SFBOFS, we assume that the model datum is equal to the North American Vertical Datum of 1988 (NAVD88) minus 0.955 m (this resultant level is close to the MSL at open ocean boundary).Therefore, an additional field, model datum minus mean sea level, was developed.In San Francisco Bay, NAVD88 data were available from Point Reyes up to the river inflow locations.
A program was developed to access the VDATUM database and to interpolate onto the SFBOFS grid the following four datum fields: MLLW to MSL, MLW to MSL, MHHW to MSL, and MHW to MSL.In addition, the specification of the model datum (MD) to MSL allows the model predicted water level results to be presented with respect to all of the tidal datums.MSL, MLLW, NAVD88 and MSL-MD of key stations are listed in Table 2.
Note that MSL-MD difference increases from the Bay entrance to Antioch and Rio Vista.The MSL at Antioch and Rio Vista are 0.20 and 0.22 m above model datum, respectfully.The digital relationships among the different tidal datums, the model datum and NAVD88 are helpful in correctly comparing model results with measured water level data.The Experiment 10 water level response at Port Chicago with respect to MLLW is shown in Figure 5.Note by using the stage boundary condition with the offsets in Experiment 10, the agreement with observations is reduced from 19 cm in Figure 4 with the flow boundary condition to 9 cm RMSE.This is due in large measure to the improvement in the simulated tidal range.

Extensive Tidal Calibration
For further calibration, the model setup used for the short term tidal experiment was used over an extended 19-month simulation from April 1979 through October 1980.Meteorological forcings were specified by setting the wind speed to zero and the sea level atmospheric pressure to 1013 mb over the entire model domain.A nudging of both salinity and temperature to specified climatological values was used along the open ocean boundary.The nineteen month simulation was completed in 38 segments of approximately 15 days' duration each, with each segment restarted from the previous segment's final fields.
In Table 3, simulation segment results for water surface elevation are compared respectively to harmonic predictions in terms of RMS error and Willmott relative error [39], which is given by <(abs(Y-X)) 2 >/<(abs(Y-<X>)+abs(X-<X>)) 2 >, with Y the model prediction and X the observation.Station locations can be found in Figure 6.In addition, model and predicted means are compared with respect to station MLLW.In general, the water level RMS errors do not exceed 15 cm and are consistent from month to month from Port Chicago in Suisun Bay through San Pablo and mid-Bay regions, as well as in the offshore and southern regions of San Francisco Bay.At Coyote Creek, at the southern end of South Bay, while the means are in close agreement, the RMS errors range from 13 to 22 cm and often exceed 15 cm.The adjustment of the bottom friction over salt marsh regions undergoing wetting and drying may need further consideration.In Table 4, principal component direction currents at mid layer (k = 10) are compared respectively to harmonic predictions in terms of RMS error and Willmott relative error.In addition, model and predicted mean currents are given.Current amplitude RMS errors are consistent from month to month and are generally less than 35 cm/s.Willmott relative errors are less than 10% except at C-33.
A more formal skill assessment has been performed in two parts.In part one, harmonic analysis was used to compare water level and principal component current strengths for the M 2 , S 2 , N 2 , O 1 , and K 1 tidal constituents.NOS accepted harmonic constants are compared with tidal simulation results in Table 5. Favorable comparisons were obtained for all constituents at all stations.In Table 6, model principal component current strengths are compared with NOS harmonic constants.Again, comparisons are favorable for both amplitude and phase at most stations except at Station C-18 for the M 2 amplitude.In part two, model and predicted means, root mean square error, standard deviation of the error, and central frequency (at reference levels of 15 cm for water level and 26 cm/s for current) were considered.In Table 7, water level skill assessment results are given with favorable comparisons exhibited for means and RMSE at all stations with the exception of Coyote Creek, where the water level error exceeded 15 cm 33.9% of the time.In Table 8, principal component current strength skill assessment results are shown with favorable results observed at most stations except again at Station C-18.The heat flux algorithm generates no excessive temperatures and produces accurate seasonal heating and cooling [22].No comparisons with observed salinity are made, since meteorological forcings are not included.However, the simulated salinity gradients are reasonable and a density front is present with the inclusion of the freshwater inflows [22].The salinity structures through the entrance are in line with climatological values.In the next section, all forcings will be turned on for the hindcast simulation and the model skill assessment will be conducted for further validation.

Hindcast Validation
An extended 19-month hindcast model validation was performed with complete meteorological forcings.The details of the model setup can be found in Section 3.During the simulation period, RMS wind speed errors are less than 5 m/s with direction RMS errors order 50 degrees.For sea level atmospheric pressure, the RMS errors are near 2 mb.For the offshore temperature and salinity, a zero gradient boundary condition is used.Since the meteorological forcings are at 3 h intervals, the effect of the sea breeze may not be completely captured.
The results are presented in 15 day increments in Table 9 for water levels.There are fewer stations available with measured data for comparison than for the tidal calibration.In addition, there are gauge datum issues at several water level stations.Generally, the water level RMS errors do not exceed 15 cm and are consistent from month to month in almost all regions.At Point Reyes, there are issues with the data, which cause errors in the subtidal water level forcings for several months indicated as blanks.
As shown in Table 10, current amplitude RMS errors are consistent throughout the period and are generally less than 35 cm/s.The salinity response is summarized in Table 11.Generally, the model salinity was in agreement with the observations at most of the stations.However, it was overestimated in the northern portion of San Pablo Bay and throughout Suisun Bay.This is believed to be due to the fact that the river subtidal water levels were not included since no measured river stage data were available.As a result, the model results could not correctly reflect the freshwater runoff during the high flow months when substantial river subtidal levels were present.This in effect, limited the amount of freshwater entering the Bay through the Delta.The temperature response is summarized in Table 12 and exhibited a normal seasonal response, but in October 1980 there was some evidence of overheating by about 2 °C in Suisun Bay.
In addition to the validation in terms of RMS errors, the NOS skill assessment criteria [40,41] are also applied to the hindcast.We show in Table 13 the results at some of the major water level stations.Additional model skill assessment results for currents, salinity, and temperature are given in [22].Generally, the skill assessment indicates that most water-level related statistical parameters pass the NOS skill assessment criteria for different scenarios, and that amplitudes and epochs of major harmonic constituents such as M 2 , S 2 , N 2 , K 2 , K 1 , O 1 , P 1 , and Q 1 from the tide-only scenario simulation are very close to the observed values at almost all stations.
Most of CF (Central Frequency), NOF (Negative Outlier Frequency), POF (Positive Outlier Frequency), MDNO (Maximum Duration of Positive Outliers), and MDPO (Maximum Duration of Positive Outliers) either pass or are close to the criteria at the Bay current stations for not only the tide-only scenario but also the hindcast scenario, since tidal current dominates the signal in San Francisco Bay region.See Schmalz [22] for more complete definitions of the skill assessment parameters.
The tidal and hindcast simulations indicate that the SFBOFS runs robustly and that the results are in acceptable agreement with the measurements.The model package was therefore loaded into the COMF-HPC on NCEP's high performance computers to perform semi-operational nowcast/forecasts.

Semi-Operational Nowcast/Forecast Simulation
The SFBOFS runs four cycles each day.In each cycle, the model performs a six hour nowcast followed by a 48 h forecast.During the model preparation process, the COMF-HPC automatically searches for and obtains the necessary observed data and other model (e.g., NAM4 and RTOFS) generated data to obtain the required forcings.For the forecast, along the open ocean boundary, water levels are specified as a superposition of the tide predictions and the subtidal water level forecast.Note the nowcast adjustments are maintained for the forecast period for water level, salinity and temperature open boundary conditions.
For both nowcast and forecast, the most recent NAM4 results in NCEP's data tank are input into COMF-HPC to get the necessary input surface forcings.
The methodology to treat the Sacramento and San Joaquin River forcings in nowcast and forecast scenarios is different because no river stage subtidal signals are available in the forecast period.Even during the nowcast, stage data are not necessarily available.COMF-HPC uses the following approach to handle this.
Subtidal river stage data adjustment is performed on the boundary nodes of the two rivers.Real-time observed stage height data from USGS 11337190 Station are taken for the San Joaquin River nodes adjustment and the data from USGS 11455420 are used for the Sacramento River nodes.The real challenge though is how to determine the subtidal water level time series for the whole nowcast and forecast time window.
As shown in Figure 7, the green curves indicate the subtidal stage height time series, which is computed as the direct water level measurement minus the tidal prediction.The vertical black time line is the current run cycle time, for example 12Z.Since the cron job is launched after the cycle time (after the NAM4 and RTOFS forcings of the same cycle are obtained), the USGS river stage reading end time, RT(end), is always on the right side of the black time line.The reading start time RT(I), however, can be on either side of the nowcast start time, ZetaT(I).
The ultimate goal is to obtain the subtidal stage height for the whole nowcast and forecast period, the time between ZetaT(I) and ZetaT(end).For the upper case in Figure 7, when RT(I) is later than ZetaT(I), we assume that stage height between ZetaT(I) and RT(I) equals the height at RT(I).For model stability, the subtidal stage height from RT(I) to RT(end) is decomposed into -mean‖ and -fluctuating‖ parts.The -mean‖ is indicated by the horizontal black line in Figure 7, and the -fluctuating‖ part is in green.The -fluctuating‖ part at RT(end) is ramped off lineally to zero in the next six hours.The -fluctuating‖ part in the rest of forecast time period is therefore taken as zero.In other words, the subtidal stage height in this period is the -mean‖ from RT(I) to RT(end).

Semi-Nowcast/Forecast Results
The SFBOFS semi-operational nowcast and forecast model assessment period started from 10 March 2013 and continued to 10 June 2013.The results from these simulations were concatenated into continuous time series for analysis using the NOS skill assessment software [17].The model ran robustly in the whole assessment period.Generally, the results of water level, current, temperature and salinity agree well with observations, and CF, NOF, POF, MDNO, MDPO, WOF and other statistical variables pass the criteria in both nowcast and forecast scenarios.Figure 8, as an example, shows the agreement of model results and observation of water level at three major stations.Refer to Peng and Zhang [42] for complete model skill assessment results at all stations for the water level, current, salinity, and water temperature.Semi-nowcast/forecast model performance is statistically shown in Figure 9.The Taylor diagrams [43] indicate that the water level results are better than the water temperature and salinity.Water level correlation coefficients at all stations are higher than 0.98, while the salinity correlation coefficient at S1 is only about 0.50 for both nowcast and forecast scenarios.The normalized modeled standard deviation at all stations is close to 1.0 for water level, but it is higher than 2.0 for some stations for salinity.Similar to the hindcast scenarios, as mentioned previously, the water level performs the best followed by water temperature and salinity.One should note that the RMSD value shown in these normalized Taylor diagrams needs to be multiplied by its corresponding measured standard deviation as listed in Tables 14-16 to get its real value.The semi-nowcast/forecast results can be found on the SFBOFS web page [44].To serve the San Francisco Bay maritime community, the SFBOFS provides users with nowcast and forecast guidance for water levels, currents, water temperature, and salinity out to 48 h, four times per day.The SFBOFS model domain on the web is divided into two separate subdomains (the San Francisco Bay and the San Francisco Bay Entrance), allowing users to focus on their area of interest.Nowcast/forecast animations of each of the two subdomains as well as time series at over 50 locations are available for winds, water levels, currents, temperature, and salinity.
Figure 10 is a snapshot from nowcast salinity animation of the larger subdomain at 0600 PST of 2 December 2013.Figure 11 illustrates that model salinity results agree well with the measurement at locations where Sacramento and San Joaquin Rivers have noticeable effect on salinity distributions.Meanwhile the available measurement at Port Chicago indicates, as shown in Figure 12, that the water level nowcast is also in good agreement with observations.The satisfying model results for both water level and water salinity near the two rivers are largely due to the fact that river stage boundary conditions have been employed.
The SFBOFS webpage offers not only the latest model output graphics as shown in Figures 10-12, but also links where users can get access and download one-year historic output files (in NetCDF format) through CO-OPS's OPenDAP and THREDDS servers.

Conclusion and Discussion
This paper details how the SFBOFS was setup, tested and extensively validated in tidal and hindcast scenarios.The performance of the model package during the three-month semi-operational nowcast and forecast using the NOS COMF-HPC is discussed.FVCOM, the core of SFBOFS, ran robustly during the trial.Amplitudes and epochs of the M 2 S 2 , N 2 , K 2 , K 1 , O 1 , P 1 , and Q 1 constituents from the tide scenario simulation are very close to the observed values at all stations.NOS skill assessment and RMS errors of all variables indicate that most statistical parameters pass the assessment criteria for both hindcast and nowcast/forecast scenarios and model outputs have good agreements with the measurement.We have to note that OTIS Regional Tide Solutions harmonic analysis results were reduced by 5% on the open boundary.Though this ad hoc treatment ensures very good water level results, more work needs to be done to understand the dynamics behind the adjustment.
Modeled water level and salinity from Martinez to Mallard Island (see Figure 10 for locations) showed strong disagreement with measurement during hindcast period when flow river boundary conditions were employed for the Sacramento and San Joaquin Rivers.The model water level results after using stage river boundary for the two rivers were greatly improved.However, salinity disagreement still existed, though in very low occurrence, during the past ten months after the semi-operational nowcast/forecast trial period.As shown in Figure 13, the model predicted salinity at Port Chicago was in agreement with the observations from October 15-October 25.However, on October 26 and 27, the model salinity predictions at Port Chicago abruptly deviated from the observations.On the 10/27/18Z cycle, the model under predicted the salinity by up to 8 PSU in the nowcast time window.A comparison of the model river forcing water surface elevation to the USGS stage data at the two rivers showed no indication that the model stage was in error.
The location of the river boundaries is still within the tidal domain and either a stage or flow boundary condition is not entirely appropriate.In effect, the boundary location is not at the head of tide and is a tidal river with flow in both directions.In the case of a stage boundary condition, no unique stage discharge relationship exists.The stage is a function of both the discharge and the offshore subtidal water level.The imposition of the stage boundary condition yields accurate water level prediction, but is problematic for salinity, since the appropriate discharge cannot always be specified.For a flow boundary condition, since the boundary is not at the head of tide, tidal wave reflections will occur and will lead to inaccurate stage predictions in the lower delta and even at Port Chicago.
In addition, the model grid cannot represent the complex water channel system in the delta region.One can compare the real delta system in Figure 1 and the model grid in Figure 2.While previous work [15] has used a 20 m deep rectangular -false delta‖ to produce the appropriate tidal prism of the unresolved area, this approach was not used in the present SFBOFS, since it was felt that the entire delta region may need to be represented as discussed by MacWilliams et al. [45].This further effort initially considered outside the scope of the SFBOFS is now being considered to improve the salinity prediction in the lower delta and to also potentially provide additional navigation guidance to the Ports of Stockton and Rio Vista.As an interim measure, a data assimilation scheme is being considered within the present SFBOFS, to correct the model salinity predictions from Martinez to Rio Vista on the Sacramento River and Antioch on the San Joaquin River at the start of each nowcast/forecast cycle.Future work will also consider the extension of the offshore boundary to include the Farallon Islands, which will allow for a more accurate specification of the offshore water level and current boundary conditions.

Figure 1 .
Figure 1.The bathymetry of San Francisco Bay Operational Forecast System and the major gauge stations.The definition of -delta‖ in this paper is the region with complicated water channels to the east of Antioch and Rio Vista.The locations of the three major bays are indicated.Note: the domain in this figure is a little larger than the model grid domain as shown in Figure 2.

Figure 2 .
Figure 2. The SFBOFS grid structure and the open boundaries.Measured river flow data from USGS are used as river forcings for the five small rivers in blue.Measured river flow or river stage data at Rio Vista and Antioch may be used as river forcings, respectively, for Sacramento River and San Joaquin River.Technically, if river stage data are employed, the grid points across these two major rivers are treated the same way as those on the Pacific Open Boundary.

Figure 3 .
Figure 3.The SFBOFS grid resolution structures are depicted.The element length sizes range from order 400-1700 m along the ocean boundary to order 200 m at the Bay Entrance as shown in the top panel.Within the Bay region the element lengths are reasonably uniform of order 150 m as shown in the lower panel with finer resolution around a few small islands of order 50 m.

Figure 4 .
Figure 4. Comparison of modeled versus predicted water level at Port Chicago with flow boundary condition over the period 1-15 April 1979.

Figure 5 .
Figure 5.Comparison of modeled versus predicted water level at Port Chicago with stage boundary condition and 5% harmonic amplitude reduction over the period 1-15 April 1979.

Figure 6 .
Figure 6.NOS and USGS Historical Circulation Survey Water Level, Current, Salinity, and Temperature Stations.Note current meters were collocated with conductivity-temperature sensors.Note the location of Point Reyes is shown in Figure 1.

6. 1 .
COMF-HPC Generated Input Forcings For the nowcast, the subtidal water levels along the open ocean boundary are determined using an adjustment of the Global RTOFS (G-RTOFS) latest hourly subtidal forecast guidance.The adjustment is determined by averaging the hourly subtidal anomalies at Point Reyes (NOAA gauge) over the previous six-hour nowcast period and ramping the forecast subtidal values to the adjustment.The astronomical tide is determined from the tidal constituent netCDF file and the application of the latest node factor and equilibrium argument values at six-minute intervals.The total open ocean boundary six-minute water level values are the sum of the adjusted subtidal levels and the predicted tidal values at each boundary grid point.Salinity and temperature along the open ocean boundary are obtained from the adjusted G-RTOFS forecast guidance.The adjustment is determined by averaging the salinity and temperature anomalies at San Francisco (NOAA gauge) over the previous six-hour nowcast period and ramping the nowcast values to the adjustment.

Figure 7 .
Figure 7. Diagram on how measured river stage height is used in COMF.

Figure 8 .
Figure 8.The comparison of modeled versus observed sea levels at three stations in April 2013.The station locations can be found in Figure 1.

Figure 9 .
Figure 9. Normalized Taylor diagrams of water level, surface temperature and surface salinity for nowcast and forecast scenarios.S1, S2, S3, etc. are station series numbers.Si of water level, temperature and salinity does not necessarily indicate the same station.The modeled standard deviation of each variable at each station is normalized by its corresponding observed standard deviation.The data are from October 2013.

Figure 11 .
Figure 11.The nowcast/forecast water surface salinity versus measurement at stations where the Sacramento and San Joaquin Rivers have noticeable effects on salinity distrubution.

Figure 12 .
Figure 12.The nowcast/forecast versus observed water levels at Port Chicago.

Figure 13 .
Figure 13.Modeled salinity strays away from observation with time.

Table 1 .
Delta Inflow Bottom Friction Experiment Summary.The scale factor was used to multiply bottom roughness in model domain above Carquinez Strait.The tapered scale factor ranges from 1 to the full value in a linear fashion from Carquinez Strait to the river inflows based on longitude.The bottom roughness sets are given in the second table.The HA amplitude reduction corresponds to reducing the amplitudes of the offshore boundary harmonic constants.

Table 2 .
Water Level Vertical Datums.Note tidal datums and NAVD88 are with respect to gage zero.Model Datum (MD) is given with respect to MSL.Note at the up estuary stations, MSL is above the model datum, while at the entrance to the Bay, MSL and the model datum are coincident.Using the table, it is possible to determine MLLW with respect to MD.

Table 3 .
Water Surface Elevation Tidal Calibration: April 1979-October 1980.For each box, the first column of values corresponds to the first 15 days of the month, with the second column of values denoting the remaining portion of the month.Within each column: Row 1 corresponds to the RMSE in cm.Row 2 corresponds to the Willmott Relative Error in percent.Row 3 corresponds to the model mean in cm with respect to MLLW with Row 4 denoting the predicted water level mean in cm relative to MLLW.

Table 4 .
Principal Current Direction Mid-Level Current Speed Tidal Calibration: April 1979-October 1980.For each box, the first column of values corresponds to the first 15 days of the month, with the second column of values denoting the remaining portion.Within each column: Row 1 corresponds to the RMSE in cm/s.Row 2 corresponds to the Willmott Relative Error in percent.Row 3 corresponds to the model mean in cm/s.Note the predicted mean current speed is zero.

Table 5 .
Harmonic Analysis of 19-month tidal simulation water level comparison to NOS accepted tidal constituents.Accepted constituent values of amplitudes (m) and phase (degrees) minus model predictions.Note negative amplitude differences denote an under prediction of tidal water level constituent amplitudes, while positive phases denote a lag in tidal water level constituent propagation.

Table 6 .
Harmonic Analysis of 19-month tidal simulation principal current direction current comparison to NOS accepted tidal current constituents.Accepted constituent values of amplitudes (m/s) and phase (degrees) minus model predictions.Note negative amplitude differences denote an under prediction of tidal current constituent amplitudes, while positive phases denote a lag in tidal current constituent propagation.Along with station id measurement, depths are given in parenthesis with observed/model principal component directions given following the depth information.

Table 7 .
Nineteen-month Tidal Simulation Water Level Skill Assessment Results.Results are with respect to MSL.RMSE is root mean square error, SD is standard deviation of the error; e.g., model minus observation, and CF is central frequency of the errors with respect to a reference level of 0.15 m.

Table 8 .
Nineteen-month Tidal Simulation Principal Component Direction Current Strength Skill Assessment Results.Note measurement depth given in parenthesis followed by observed/model principal current direction.RMSE is root mean square error, SD is standard deviation of the error; e.g., model minus observation, and CF is central frequency of the errors with respect to a reference level of 0.26 m/s.

Table 9 .
Water Surface Elevation Hindcast Validation: April 1979-October 1980.For each box, the first column of values corresponds to the first 15 days of the month, with the second column denoting the remaining portion of the month.Within each column: Row 1 corresponds to the RMSE in cm.Row 2 corresponds to the Willmott Relative Error in percent.Row 3 corresponds to the model mean in cm relative to MLLW.Row 4 denoting the observed water level mean in cm with respect to MLLW.

Table 10 .
Current Speed Hindcast Validation: April 1979-October 1980.In each box, the first column of values corresponds to the first 15 days of the month, with the second column of values denoting the remaining portion.For each column: Row 1 corresponds to the RMSE in cm/s.Row 2 corresponds to the Willmott Relative Error in percent.Row 3 corresponds to the model mean in cm/s with Row 4 denoting the observed mean current speed in cm/s.Note n/a denotes not available due to incomplete observations.Along with the station id, the measurement location (in distance above the bottom in meters) is given in parenthesis.

Table 11 .
Salinity Hindcast Validation: April 1979-October 1980.In each row box, the first column of values entry corresponds to the first 15 days of the month, with the second column denoting the remaining portion.Within each column: Row 1 corresponds to the RMSE in PSU.Row 2 corresponds to the Willmott Relative Error in percent.Row 3 corresponds to the model mean in PSU with Row 4 denoting the observed salinity mean in PSU.Note n/a denotes not available due to incomplete observations.Along with the station id, the measurement location (in distance above the bottom in meters) is given in parenthesis.

Table 12 .
Temperature Hindcast Validation: April 1979-October 1980.In each box month, the first column values correspond to the first 15 days of the month, with the second column of values denoting the remaining portion.Within each column: Row 1 corresponds to the RMSE in °C.Row 2 corresponds to the Willmott Relative Error in percent.Row 3 corresponds to the model mean in °C with Row 4 denoting the observed temperature mean in °C.Note n/a denotes not available due to incomplete observations.Along with the station id, the measurement location (in distance above the bottom in meters) is given in parenthesis.

Table 13 .
Nineteen-month Hindcast Water Level Skill Assessment Results.Results are with respect to MLLWL.RMSE is root mean square error, SD is standard deviation of the error; e.g., model minus observation, and CF is central frequency of the errors with respect to a reference level of 0.15 m.

Table 14 .
Observed water level standard deviations (m) at selected stations in nowcast and forecast scenarios.

Table 15 .
Observed temperature standard deviations (°C) at selected stations in nowcast and forecast scenarios.

Table 16 .
Observed salinity standard deviations (PSU) at selected stations in nowcast and forecast scenarios.