OCS Study 

BOEM 2015-035

Data Assimilative Hindcast 

for the Gulf of Mexico

US Department of the Interior 

Bureau of Ocean Energy Management

Headquarters, Sterling, VA
i 

OCS Study 

BOEM 2015-035

Data Assimilative Hindcast 

for the Gulf of Mexico

Eric P. Chassignet

Center for Ocean-Atmospheric Prediction Studies, Florida State University, 

Tallahassee, FL

Ashwanth Srinivasan

Tendral LLC, Miami, FL

Prepared under BOEM Award M12AC00019

US Department of the Interior 

Bureau of Ocean Energy Management 

Headquarters 
ii

DISCLAIMER

Study concept, oversight, and funding were provided by the US Department of the Interior, 

Bureau of Ocean Energy Management, Environmental Studies Program, Washington, DC, 

under Contract Number M12AC00019. This report has been technically reviewed by BOEM 

and it has been approved for publication. The views and conclusions contained in this document 

are those of the authors and should not be interpreted as representing the opinions or policies of 

the US Government, nor does mention of trade names or commercial products constitute 

endorsement or recommendation for use. 

REPORT AVAILABILITY 

To download a PDF file of this Environmental Studies Program report, go to the US Department 

of the Interior, Bureau of Ocean Energy Management, Environmental Studies Program 

Information System website and search on OCS Study BOEM 2015-035. 

This report can be viewed at select Federal Depository Libraries. It can also be obtained from the 

National Technical Information Service; the contact information is below. 

US Department of Commerce 

National Technical Information Service 

5301 Shawnee Rd. 

Springfield, VA 22312 

Phone: (703) 605-6000, 1(800)553-6847 

Fax: (703) 605-6900 

Website: http://www.ntis.gov/

CITATION

Chassignet E.P., and A. Srinivasan, 2015. Data Assimilative Hindcast for the Gulf of Mexico. 

US Dept. of the Interior, Bureau of Ocean Energy Management, Headquarters, Sterling, VA. 

OCS Study BOEM 2015-035. 46 pp. 

Cover Page Image Information: A sample snapshot of the surface temperature (o

C) from the 

hindcast for January 2, 2003. 

ACKNOWLEDGMENT

This study was funded by the US Department of the Interior, Bureau of Ocean Energy 

Management under the cooperative agreement M12AC00019. Zhen Li is the Project Officer of 

this study. 
iii

Data Assimilative Hindcast for the Gulf of Mexico

Contents

Contents iii

List of Figures iv

List of Tables iv

List of Acronyms v

Executive Summary vi

1. Introduction 1

2. Methodology 2

2.1 The Ocean Model and its Configuration for the Gulf of Mexico. 2

2.2 Observations for Providing Data Constraints 2

2.3 The Ocean Data Assimilation System 3

3. Initial Conditions for the Hindcast 5

4. Qualitative Evaluation of the 2003-2012 Hindcast 6

5. Quantitative Comparisons 9

5.1 Comparisons with the Along Track Altimeter Data 9

5.2 Comparison with SST Data 11

5.3 Comparison with In-Situ Profiles 12

5.4 Comparison with Surface Drifting Buoys 15

5.5 Verification with Oil Spill Trajectories 17

6. Discussion and Conclusions 17

7. Project Deliverables 18

References 20

Appendix A. Optimizing the Atmospheric Forcing 22

A-1. Radiation Flux Corrections for Thermal Forcing 22

A-2. Wind Speed Correction Based on Linear Regression to QuickSCAT Winds 23

Appendix B. Observation Error Specification and processing 25

Appendix C: Tendral-Statistical Interpolation system (T-SIS) implementation details for HYCOM 27

Appendix D: System Documentation 29

D-1. Overview 29

D-2. Contents of the Package 29

D-3. Compiling the package 29

D-3.1 Required Compilers and Libraries 29
iv

D-3.2 Compiling the HYCOM ocean model 30

D-3.3 Compiling the T-SIS Assimilation Codes 31

D-4. Setting up the 2010 Hindcast 31

D-5. Reference Model Inputs Parameters 33

D-5.1 HYCOM Reference blkdat.input 33

D-5.2 T-SIS Reference Namelist 36

 

List of Figures

Figure 1. The multiscale assimilation methodology. 5

Figure 2.Initial model state on January 1, 2003 derived from the Global HYCOM Reanalysis. 6

Figure 3. Loop current state on July 9, 2003 (sea surface height). 7

Figure 4. Loop current state on July 9, 2003 (ocean color and ocean current speed). 8

Figure 5. Time series of rms errors between model hindcast and assimilated sla observations. 10

Figure 6. Altimetry and model based estimates of the mean and variance of the sea surface elevation. 11

Figure 7. Sea surface tempreature (°C ) for july 9, 2003. 12

Figure 8.Conductivity temperature depth (CTD) stations from the wod 2009 for the 2003-2012 time frame. 13

Figure 9. Mean bias and rmsd profiles for temperature and salinity. 14

Figure 10. Location of the P3 airborne expendable bathythermograph (AXBT) observations collected during the 

DwH incident. 14

Figure 11. Temperature bias (°C ) and RMSD at AXBTS’ location. 15

Figure 12. Location of velocity estimates from the drifters. 16

Figure 13. Velocity RMSD and observed standard deviation; model drifter velocity correlation. 16

Figure 14. Location of the particles one month after release at the wellhead. 18

Figure 15. Scatter plot between monthly mean QSCAT winds and the NARR winds. 24

Figure 16. Schematic of the statistical interpolation algorithm used in the T-SIS assimilation scheme. 28

List of Tables

Table 1. Overview of the available observations in the Gulf of Mexico 3

Table 2. Timing of separation/reattachment and shedding events as reproduced by the hindcast and CCAR 

altimetry 8

Table 3. Total number of profiles available in the World Ocean Database for the 2003-2012 hindcast 13

Table 4. Summary of the assumed standard deviations of instruments errors. 25
v

List of Acronyms

AOML Atlantic Oceanographic and Meteorological Laboratory

AVHRR Advanced Very High Resolution Radiometer

CCAR Colorado Center for Astrodynamics Research

CTD Conductivity, Temperature, and Depth

DwH Deepwater Horizon

ECMWF European Center For Medium Range Weather Forecasting

ERM Exact Repeat Mission

ERS European Remote Sensing

GDEM Global Digital Elevation Map

GHRSST Group for High Resolution Sea Surface Temperature

GoM Gulf of Mexico

GFO Geosat Follow-on

HYCOM HYbrid Coordinate Ocean Model

MODIS Moderate Resolution Imaging Spectroradiometer

NARR North American Regional Reanalysis

NCEP National Centers For Environmental Prediction

NAVOCEANO Naval Oceanographic Office

NOAA National Oceanic & Atmospheric Administration

NODC National Oceanographic Data Center

OSTM Ocean Surface Topography Mission

QuickSCAT Quick Scatterometer 

RMS Root Mean Square

RMSD Root Mean Square Deviation

SLA Sea Level Anomaly

SSH Sea-surface height

SST Sea Surface Temperature

T/P Topex/Poseidon

T-SIS Tendral–Statistical Interpolation System 

T/S Temperature/Salinity
vi

Executive Summary

Accurate estimates of ocean currents are critical for oil spill risk analysis. As the government 

agency responsible for oil spill risk analysis, the Bureau of Ocean Energy Management (BOEM)

requires estimates of ocean currents for oil spill risk analysis in the Gulf of Mexico. The 

cooperative agreement, M12AC00019, was initiated to produce a validated 2003-2012 Gulf of 

Mexico ocean state hindcast. Toward this end, we have used the latest version of the HYbrid 

Coordinate Ocean Model (HYCOM) and an ensemble based data assimilative framework to 

estimate the ocean state variables necessary for oil spill risk analysis. Along-track sea level 

anomalies, gridded sea surface temperature, in-situ and climatological profiles of salinity and 

temperature were assimilated to provide data constraints to the ocean model. An assessment of 

the hindcast shows robust model skill in reproducing the various features of the circulation in the 

Gulf of Mexico. As expected for a model assimilating sea level anomalies, there is a good 

agreement between the hindcast and independent estimates of Loop Current behavior inferred 

from altimetry and ocean color both in terms of temporal and spatial characteristics. 

Furthermore, the assimilation of a limited number of in-situ and climatological profiles of 

temperature and salinity is effective in controlling bias and reduces errors in the vertical density 

structure in the model. Comparison with both assimilated and unassimilated in-situ observations 

indicate that the upper ocean temperature and salinity errors are within 1.5°C and 0.3 psu for the 

hindcast time period. Comparisons with independent drifter derived velocity data indicate 

significant correlation (>0.5) between the drifter data and model velocities. Simplified oil spill 

simulations of the Deep Water Horizon event captures the main features of surface plume 

indirectly verifying the underlying model circulation fields. Overall, the hindcast dataset is 

consistent with what is observed during the 2003-2012 time frame and is appropriate for a range 

of applications including oil spill risk analysis. This report details the methodology and data 

used for the hindcast and includes an assessment of the 2003-2012 hindcast. The complete model 

hindcast dataset, model code, documentation, and a test case for the year 2010 is provided to 

BOEM along with this report. 
1

1. Introduction

The Gulf of Mexico (GoM) region accounts for 30% of the total oil production in the US. At 

present, there are more than 3,500 oil platforms generating more than 1.7 million barrels of oil 

per day. Oil production in the GoM is expected to continue at these levels into the next decade 

and thus carries the risk of potential oil spills into the future. Oil spill risk analysis is critical for 

evaluating the potential environmental impacts associated with the oil and gas development in 

the Gulf and underpins all preparations for spill response. In these analyses, ocean currents and 

their statistics are used for simulating trajectories to infer potential pathways and rates of 

contamination in case of a spill. A necessary requirement for a robust and scientifically rigorous 

risk assessment is an accurate knowledge of the statistics of powerful, complex and highly 

variable ocean currents. One way of obtaining such estimates is by means of data assimilative 

hindcasts. Such an hindcast blends observations with a dynamical circulation model to provide 

regular three dimensional estimates of ocean currents taking into account the errors in both 

observations and models.

Over the years, ocean models used for studying the GoM circulation have become increasingly 

sophisticated (Oey et al., 2005). However, most models have been configured for process studies 

and do not address practical needs adequately. For many practical applications, a long time time 

series of model outputs from realistically configured ocean models, forced with observed or 

analyzed winds and assimilating observations to incorporate data constraints is desirable. Data 

constraints placed by the assimilation process keeps the evolving model state close to the 

observed state during the course of the hindcast. A few models such as the US Navy’s GoMHYCOM, used for routine operational forecasting, satisfy many of these requirements. However, 

improvements and changes in the operational system over time means that the time series of the 

outputs from such a system is not consistent temporally and/or spatially and might be of uneven

quality. Using operational outputs for applications such as oil spill risk analysis is generally not a 

straightforward task. A reanalysis or a long-term model run with a single model configuration

and data assimilation method can, however, provide a consistent dataset for a wide range of 

applications. Kantha et al. (2005) describes such an analysis and point out the skill and utility of 

the hindcasts in complementing the observations. For their analysis, only remotely sensed 

altimeter data was used but vertical profiles were not used and no constraints were placed on the 

vertical density structure. With the risk of undersea oil spill, dramatically illustrated by the 

Deepwater Horizon (DwH) event, it is important to accurately estimate both surface and 

subsurface currents as well as the vertical density structure. This is especially important for oil 

spill risk analysis since the subsurface aspects of the circulation controls the transport and fate of 

oil trapped in the water column. Assimilation of in-situ profiles of salinity and temperature is 

therefore required to reduce model errors in the water column. In situ data in the Gulf are sparse 

and unevenly sampled in both space and time, but, even so, assimilating a limited number of 

vertical profiles of temperature and salinity along with remotely sensed altimeter data is expected 

to improve the modeled circulation. For this study, we have set up a 10-year hindcast (2003-

2012) that assimilates along-track sea level anomalies, sea surface temperature ,and profiles of 

salinity and temperature. In the following sections of this report, we first describe the ocean 

model configuration, the observations, and the assimilation methodology used for the hindcast. 

We then provide an assessment of the 2003-2012 hindcast. 
2

2. Methodology

2.1 The Ocean Model and its Configuration for the Gulf of Mexico.

We use the HYbrid Coordinate Ocean Model (HYCOM, http://hycom.org), a model that is 

widely used by the oceanographic community and and is skillful at representing physical

processes in the Gulf of Mexico. The centerpiece of HYCOM is a generalized vertical coordinate 

designed to quasi-optimally resolve vertical structure throughout the ocean (Bleck, 2002; 

Chassignet et al., 2006). Typically, it is isopycnic in the stratified interior to minimize spurious 

diapycnal mixing, level or pressure coordinates near the surface to provide resolution in the 

surface mixed layer, and terrain-following in coastal regions to accurately represents flowtopography interactions. The nature of the model’s vertical layers vary both in space and time 

with the background state at any given instant and is defined by a vertical grid generator. The 

coordinate system adopted in HYCOM allows for the use of advanced subgrid scale 

parameterizations in the ocean mixed layer while retaining the advantages of an isopycnal model 

in the ocean interior (Bleck, 2002, Chassignet et al., 2003, Chassignet et al., 2006). HYCOM is 

maintained as a portable code that runs efficiently on most modern architectures. The code can 

be used on shared memory, distributed memory and hybrid computer systems. The parallel 

model outputs are reproducible and can be saved in many formats including NetCDF. 

We have adopted a model configuration based on a 1/25o GoM-HYCOM system that is used in 

near real time by the US Navy and by the scientific community for process studies. The model 

bathymetry is derived from a recent high resolution bathymetry available from Florida State 

University (see report by Velissariou (2011) for details). For the surface forcing fields, we use 

the North American Regional Reanalysis (NARR) atmospheric product (Messinger et al., 2006). 

Surface fluxes of heat and momentum are calculated through bulk formula using fields of wind 

stress, 10 m wind speed, 2 m air temperature, 2 m water vapor mixing ratio, precipitation, and 

incoming solar radiation. A correction to the NARR derived radiation fluxes were derived by 

comparing model surface temperature with Advanced Very High Resolution Radiometer 

(AVHRR) estimates of surface temperature fields. The NARR wind speed was corrected by 

regression to QuikSCAT wind speeds. More details on these corrections are listed in Appendix 

A. Monthly climatological Gulf of Mexico river data from USGS database 

(http://waterdata.usgs.gov/nwis/sw) and the RIVDIS database (http://www.RivDis.sr.unh.edu) 

were specified for the freshwater discharge (treated as a virtual salt flux). The hindcast 

simulation uses bi-weekly climatological lateral open boundary conditions derived from a North 

Atlantic climatological run. 

2.2 Observations for Providing Data Constraints

Observations are at the heart of the data assimilative system since they provide the necessary 

data constraints crucial for error reduction and validation. We have used several sources of 

observations, including remote sensed satellite Sea Surface Temperature (SST) and Sea Level 

Anomalies (SLA), in-situ temperature and salinity profiles obtained from ship surveys, 

moorings, profiling floats, and gliders. The characteristics of the data and their availability are 

listed in Table 1 below. 
3

Table 1. Overview of the available observations in the Gulf of Mexico

Type Provider/Source Frequency Spatial 

Characteristics

1 SLA Collecte Localisation Satellites,

www.aviso.oceanobs.com

Daily Along track

2 Foundation 

Temperature

ftp://podaac-ftp.jpl.nasa.gov/allData/GHRSST Daily Gridded

3 0.25° SST analysis http://www.ncdc.noaa.gov/sst/ Daily Gridded

4 In-situ 

Temperature/Salinity

NOAA/NODC World Ocean Data Base 

www.nodc.noaa.gov/OC5/WOD09/pr_wod09.html

Irregular Point

5 In-situ 

Temperature/Salinity

NOAA/NODC Deepwater Horizon response data

http://www.nodc.noaa.gov/General/DeepwaterHori

zon/support.html

Irregular Point

Remotely sensed sea level anomalies (SLA) and SST as well as in-situ temperature/salinity (T/S)

profiles are considered to be the most reliable and are systematically assimilated. Along-track 

SLA are available for the hindcast period from multiple operational satellite altimeters,

Topex/Poseidon (T/P), Jason 1,2, Envisat, Geosat Follow-on (GFO) and Cryosat. The data are 

available from Collecte Localisation Satellites (CLS) from January 1993 to present in a delayed 

time quality controlled mode and is the primary data set for constraining the model.

For the sea surface temperature, we use the 0.25° daily SST analysis available from NOAA. Insitu profiles of temperature and salinity from the World Ocean Atlas 2009 were assimilated. Data 

from cruises during the DwH incident were also used to constrain the sub-surface density 

structure. The other observations – acoustic Doppler current profilers (ADCPs), high frequency 

radar and drifters – were not assimilated. Among these, we use only the drifter data for 

validations since data from the other two platforms are of unknown quality. Although all 

observations are obtained from delayed mode quality controlled datasets, additional checks and 

careful preprocessing of the data was undertaken to ensure that there are no obvious erroneous 

values in the data stream since such values can lead to the persistence and propagation of these 

errors in a data sparse environment. Further details on quality control procedures and observation 

error variances are listed in Appendix B.

2.3 The Ocean Data Assimilation System

We use the Tendral–Statistical Interpolation System (T-SIS), a statistical interpolation package 

optimized for the Lagrangian vertical coordinate of HYCOM (Appendix C). An application of 

this package for the Gulf of Mexico is documented in Halliwell et al. (2014).

The assimilation scheme is based on multivariate linear statistical estimation wherein the best 

linear unbiased estimate of the state of the ocean (the analysis ) is obtained by updating the 

previous model forecast using

 (1)

a x

f x

( ), a f f x x K y Hx =+ −
4

where represents the observations to be assimilated, is the observation operator, is a 

matrix of optimization parameters often called the gain matrix. The Gauss-Markov formula 

prescribes a gain matrix that is optimal in a least-square sense (e.g. Bennett, 1992; Wunsch, 

1996):

 (2)

where is the forecast error covariance matrix, is the observation error covariance matrix, 

and superscript T denotes matrix transpose. Formally, is the covariance matrix of 

the forecast error given , where E represents an ensemble average and is the 

true state of the ocean. The forecast is assumed to be statistically unbiased with .

Due to the lack of accurate information on the true oceanic state, is a difficult quantity to 

determine. Numerous approximations of have been used to represent the multivariate and 

spatial correlations as accurately as possible in a numerically efficient fashion. The approach 

taken here is to prescribe the error covariance using an ensemble of model states sampled at 

different times:

 (3)

where is the m

th sample of the forecast ensemble, is the ensemble mean, and M is the 

number of samples. The underlying assumption is that the time variability can be related to error 

covariance. For the hindcast, this matrix is estimated from model states sampled from a long 

unconstrained run of the forecast model. For a given analysis day the model states are chosen 

from a 90 running window to account for the time variability. However, the magnitude of the 

true error covariance is likely to be smaller than the time variability. Therefore, the prescribed 

error covariance is scaled to realistic levels. Measurements are assumed to be uncorrelated and a 

diagonal observation error covariance matrix is used. Srinivasan et al. (2011) compared this 

approach with three other methods used operationally and found it to be well suited for 

constraining the mesoscale circulation in the Gulf of Mexico.

One novel feature of the assimilation scheme implemented here is the use of climatological 

profiles of temperature and salinity to supplement the sparse insitu data. The irregular sampling 

and the relatively small number of assimilated in-situ profiles makes it difficult to adequately 

constrain the vertical density structure with these profiles alone. This can lead to the slow growth 

of errors in the mean model climatology. Although these errors are likely to be within the normal 

climatological variances for 10 year integrations, they are best controlled before problems 

manifest. Direct assimilation of climatological profiles into eddy resolving models is known to 

suppress model variability and is therefore not desirable. Therefore ,we use a multiscale 

procedure designed to allow the assimilation of climatological profiles into eddy resolving 

models without suppressing the model variability. The method uses two sets of covariance: one 

for climatological scales and the other for mesocale. The assimilation is done in two stages, the 

y H K

( ) 1

, fT fT −K P H HP H R = +

f P R

( ) f f fT P ee = E

f f true exx = − true x

( ) 0 f E e =

f P

f P

( )( ) 1

1

( 1) ,

M T

f ffff

m m

m

M −

−

P =− − − ∑ xxxx

f

m x f x
5

climatological scale is first corrected and then the evolving mesoscale field is corrected as 

illustrated in Figure 1 below.

Figure 1. The multiscale assimilation methodology.

The climatological scale is first corrected with climatological profiles to 

suppress drift and bias before correcting the mesoscale.

3. Initial Conditions for the Hindcast

We seek an initial state that is dynamically balanced and representative of the conditions at the 

start date. Further, the initial state will ideally be obtained from a model whose mean vertical 

density structure is close to observed climatology. We derive the initial conditions for the 

hindcast by interpolating the 1/12° global HYCOM reanalysis (www.hycom.org/dataserver) to 

the 1/25° GoM domain. The surface speed and mixed layer temperature (0.2° C change from the 

surface) from the initial state on Jan 1, 2003 are shown below with the Loop Current extending 

well into the Gulf at this time. The three-year mean (2003-2005) vertical profiles of T/S from 

the global reanalysis climatology shows small biases between 300-700 m, but it is within the 

climatological variances at these depths. Overall, the upper ocean vertical density structure is in 

reasonable agreement with the Global Digital Elevation Map (GDEM) climatology (Carnes, 

2009) as illustrated below in Figure 2. 
6

 

 

Figure 2.Initial model state on January 1, 2003 derived from the global 

HYCOM reanalysis.

The upper panels show the surface speed and mixed layer temperature (

o

C) on 

January 1, 2003. The bottom panels show mean vertical profiles of temperature 

and salinity from the reanalysis when compared to the GDEM climatology. 

4. Qualitative Evaluation of the 2003-2012 Hindcast

The Loop Current is the principal dynamical feature of the GoM. It exhibits a complex behavior 

and is known to shed large anti-cyclonic eddies at quasi-regular intervals of 3 to 17 months 

which impacts almost all aspects of the circulation in the entire GoM. Therefore we begin the 

evaluation of the hindcast by a qualitative examination of the Loop Current and the associated 

circulation. There were several Loop Current eddy shedding, separation and re-attachment 

during the 2003-2012 time period. Their characteristics and timings as reproduced in the hindcast 

can be compared with the historical maps of optimally interpolated surface elevation provided by 

Colorado Center for Astrodynamics Research (CCAR). The historical CCAR product is a 

merged solution of Geosat Exact Repeat Mission (ERM), European Remote Sensing (ERS)-1, 

T/P, ERS-2, GFO, Envisat, Jason-1, Jason-2/ Ocean Surface Topography Mission (OSTM) and 

Cryosat-2 data products (Leben et al., 2002). The data includes +/-10 days of T/P, Jason-1 and 

Jason-2/OSTM data centered about the date of interest along with +/-17 days of Geosat, ERS-1, 
7

ERS-2, GFO, Envisat and Cryosat-2 data centered about the date of interest. While comparing 

the CCAR product with the hindcast, it is important to note that the model was run in a filtering 

mode and only used past data within a two-day window. We list the eddy shedding/separation 

and re-attachment events in Table 2 along with altimeter based estimates and illustrate one such 

event below. As can be seen from Table 2, the agreement is quite close between the two 

estimates. 

In Figure 3 below, we show a snapshot on such event on July 9, 2003 when an eddy is soon to be 

shed. We first compare the modeled surface elevation field for this day to the CCAR optimally 

interpolated sea surface elevation for the same day. The model derived surface elevation fields

are in good agreement with the surface elevation obtained from the CCAR optimal interpolation 

of altimeter data. The shape and extent of the Loop current extension and the cyclonic circulation 

features around the periphery of the Loop Current and in the western GoM seen in the altimetry 

image are all well represented in the hindcast. 

 

Figure 3. Loop Current state on July 9, 2003 (sea surface height).

Surface elevation from altimetry (left) and the hindcast (right) illustrating the 

qualitative similarity of the fields. Note a model mean has been added to the 

CCAR optimally interpolated sea leval anomalies to obtain the surface elevation 

field. 

While this agreement between the altimetry and model is necessary, it is not sufficient for 

assessing the realism of the model since the model assimilates the same altimeter observations 

used in generating the altimetry estimate. Chassignet et al. (2005) point out the utility of ocean 

color images as an independent means of assessing the performance of ocean models. Typically, 

the dark areas in the ocean color image represent the low-chlorophyll Loop Current and 

Caribbean water masses. The bright colors areas with higher chlorophyll content are generated 

along the coasts and are sometimes advected into the cyclonic circulation features. 
8

 

Figure 4. Loop Current state on July 9, 2003 (ocean color and ocean current 

speed).

The figure should be compared with Figure 3. The complex circulation pattern in 

the Eastern flank of the Loop Current where a cyclonic intrusion is visible in the 

color image is well represented in the surface speed image. 

The ocean color image obtained by composing images for a 10-day period centered on July 9,

2003 offers a complimentary view of the Loop Current state during this time period. The overall 

shape and the northward extent of the dark colored regions near the Loop Current region in the 

ocean color image are quite similar to the shape of the Loop current reproduced in the hindcast. 

The cyclonic circulation feature at the eastern flank is in the process of cleaving a Loop Current 

Eddy and is also well reproduced in the hindcast. The model surface speeds show the complex 

circulation in the eastern flank of the Loop Current. In general, such qualitative agreement 

between the model and observations is typical for the entire duration of the hindcast.

The frequency of eddy shedding, the northwestern penetration into the Gulf, the peeling away of 

the eddy from the Loop current due to westward propagation, the role of cyclones both small and 

large in necking down and pinching the Loop Current Eddy, the sizes of the eddies shed, their 

southwest propagation are identified as generic components of the Loop Current Eddy dynamics.

The animations of surface elevation anomalies of the CCAR fields and of the hindcast (available 

at http://stargazer.coaps.fsu.edu/boem_gom_hindcast) clearly display these components.

Table 2. Timing of separation/reattachment and shedding events as reproduced by the 

hindcast and CCAR altimetry

Overall, there is a good agreement and the shedding events agree with each others within 10 

days.

Year Loop Current Eddy Events CCAR-Altimeter Hindcast

2003 1) small eddy separation

2) shedding-reattachment

3) shedding

Day 220

Day 270-340

Day 360

Day 224

Day 275-340

Day 365

2004 1) shedding Day 240 Day 240
9

Year Loop Current Eddy Events CCAR-Altimeter Hindcast

2005 1) separation-reattachment

2) shedding

Day 63-102

Day 260

Day 55-100

Day 260

2006 1) small eddy shedding

2) shedding

Days 42,68

Day 268

Days 43,70

Day 270

2007 1) separation/reattachment

2)shedding

3) shedding

Days 105-145

Day 162

Day 318

Days 106-150

Day 160

Day 320

2008 1) small eddy shedding

2) shedding

Day 70

Day 185

Day 70

Day 192

2009 1) shedding

2) shedding

Day 63

Day 245

Day 65

Day 245

2010 1) separation/reattachment Days160-240 Days160-230

2011 1) shedding Day 230 Day 225

2012 1) shedding

2) shedding

Day 41

Day 200

Day 30

Day 180-200

5. Quantitative Comparisons

5.1 Comparisons with the Along Track Altimeter Data

The GOM hindcast assimilates daily along track altimeter data as the principal data constraint for 

the time period 2003-2012. We first compare time series of errors in sea level anomalies with the 

along track along track anomalies. Since these data were assimilated, they cannot be used for 

validation but a good fit and stable error levels are necessary to verify that the assimilation was 

successful in constraining the model. It is to be noted that the data were further processed by 

forming super-observations (multiple observations within 0.2° of each other were combined into 

one observation) before assimilation, the comparisons shown here however are with the raw 

data. Figure 5 shows the Root Mean Square (RMS) errors for six out of the 10 years during 

2003-2012 period. An average of 400 observations are assimilated daily. The mean errors are 

about 7 cm during the entire hindcast, but the errors are a bit higher during the 2008-2012 time 

frame when there are only three altimeters. Best fit is obtained when 5 altimeters are available. 
10

Figure 5. Time series of RMS errors between model hindcast and assimilated 

SLA observations.

On average, the errors are around 7 cm for the entire duration of the hindcast. 

However errors are marginally higher during periods coinciding with reduced 

number of altimeters.

Maps of the mean and the variance of the surface height fields derived from the objectively 

analyzed altimeter data (Leben et al., 2002) and the hindcast are shown below (Figure 6). Both 

maps below are generated using the daily data for the 2003-2012 time period. They show broadly 

similar features with the variability principally due to the Loop Current intrusion and westward 

propagation of the Loop Current Eddy. The mean extent of northward penetration and the 

latitude range of the southwestward propagation are quite similar between the two analyses and 

is in good agreement with maps of the Mean Dynamic Topography and the estimates of surface 

height variance provided by AVISO. 
11

Figure 6. Altimetry and model based estimates of the mean and variance of 

the sea surface elevation.

The maps were generated based on daily data. Note that the altimetry based 

estimates are generated at 0.25° resolution.

5.2 Comparison with SST Data

Our original intention was to assimilate the high resolution GHRSST analysis. Although some of 

the GHRSST products are of a higher resolution than the model, they still are very smooth and 

do not represent the fine scale features seen in the model. For example, the 1 km GHRSST image 

and the model output for July 9, 2003 are shown below (Figure 7). The Loop Current and the 

intrusion of colder waters from the shelf are clearly seen in the model, but only faintly in the 

GHRSST image. In light of this scale mismatch, it was decided to assimilate the 0.25° AVHRR 

data to correct only the large scales in the model. A smoothing operator was used on the model 

forecast to calculate innovations and data was directly assimilated at each grid point with 

covariance implicitly available from the gridded field. The AVHRR data minus model Root 

Mean Square Deviation (RMSD) time series (2003-20012) is shown in the bottom left panel of 

Figure 7.
12

 

 

Figure 7. Sea Surface Temperature (°C ) for July 9, 2003.

Multiscale Ultrahigh Resolution (MUR) SST analysis and model mixed layer 

temperature for July 9, 2003 (top panel). Time series of SST RMSD and mean 

winter months error map (lower panel). 

It can be seen that there is a seasonal cycle in the error pattern. The RMSD is higher during the

winter months and about 0.8°C while in summer it is typically about 0.5°C. A closer look reveals 

that most of the problems are in the shelf regions where there are mean errors of 2°C. 

5.3 Comparison with In-Situ Profiles

The hindcast assimilates profiles of temperature and salinity from the World Ocean Database

(WOD). There are about 500 profiles per year excluding the year 2010 when there is a lot more 

due to the intense observations during the DwH incident (Table 3). Figure 8 shows the 

distribution of the CTD station available from the WOD. The sampling is highly irregular. The 

northern shelf areas are well sampled while there are only a few stations in other areas. 

Furthermore, many of the stations in the northern regions are very shallow and too close to each 

others. We used a conservative, but consistent, thinning procedure (see Appendix B) for 

selecting the profiles to assimilate from the WOD database. 
13

Table 3. Total number of profiles available in the World Ocean Database 

for the 2003-2012 hindcast

Year Total Number of Profiles

2003 403

2004 591

2005 602

2006 431

2007 517

2008 267

2009 60

2010 2600

2011 304

2012 450

Figure 8.Conductivity Temperature Depth (CTD) stations from the 

WOD 2009 for the 2003-2012 time frame.

We compare the modeled and observed profiles and show mean bias and the RMSD over the 

water column depth in Figure 9. These statistics are calculated using all available observations 

(many of which were not assimilated). The temperature and salinity profiles indicate the 

presence of a relatively small mean bias for temperature and salinity in the upper 1000 m. The 

corresponding RMSD for temperature and salinity is generally less than 1.5°C and 0.5

respectively. In spite of the presence of these errors, we note that the mean square error is less 

than the climatological variance suggesting that the signal to noise ratio is greater than one.
14

Figure 9. Mean bias and RMSD profiles for temperature and salinity.

The analysis shown is for the year 2010 and uses all available observations (many 

of which were not assimilated). 

Next we compare the hindcast with the Airborne EXpendable BathyThermograph (AXBT) 

observations collected during the DwH incident. 

Figure 10. Location of the P3 Airborne EXpendable 

BathyThermograph (AXBT) observations collected during the DwH 

incident. 

The flights were conducted by NOAA P3 aircraft. The AXBT’s were 

supplied by BOEM, originally Minerals Management Service, U.S. 

Department of the Interior.

An indication of the vertical distribution and the time evolution of the errors is shown in Figure 

11. The error levels are similar to the ones seen in the profile above and are stable during the 17 

flight days between May 8 and July 25, 2010. Also shown are the error statistics from NRL’s 

global HYCOM reanalysis and GOM-HYCOM which also assimilated the AXBTs. The errors 

are of the same order of magnitude.
15

 

Figure 11. Temperature bias (°C ) and RMSD at AXBTs’ location.

The bias and RMSD over the upper 500 m during the 17 flight days in from May 

8 to July 25, 2010 (top panel) and time series of temperature errors are shown in 

the bottom panel.The location of the P3 observations are provided in Figure 10.

5.4 Comparison with Surface Drifting Buoys

Velocities from surface drifting buoys are not assimilated in the hindcast and thus constitute an 

independent dataset for skill evaluation and verification. We use the 6-hourly krigged velocities 

representative of 12 m depth (www.aoml.noaa.gov/phod/dac/dacdata.php) provided by NOAAAOML. A sample of the locations at which velocities are available for 2010 is shown in Figure 
16

12. The 12 m model velocities were interpolated at the drifter locations for 2007, 2010 and 

2012. 

Figure 12. Location of velocity estimates from the drifters.

There are close to 17,000 degrees of freedom. 

Figure 13. Velocity RMSD and observed standard deviation; model drifter 

velocity correlation.

Standard deviation of the drifter velocity estimates and RMS errors in velocities are 

shown in the left panel. The correlations between the drifter velocities and model 

velocities are shown in the right panel.

The RMS errors between the drifter velocities and model are most of the time smaller than the 

standard deviation of the drifter velocities indicating useful model skill. The correlation 

coefficients between the two velocity estimates are generally about 0.6 with the northward 

velocity showing slightly higher correlations. There are several thousand velocity estimates for 

each year, so these correlations are robust. A more thorough analysis should take into account 
17

the errors in direction in addition to the errors in magnitude shown above. Another issue is that 

most of the drifter locations are in the Loop Current region and it is not clear how representative 

these statistics might be in other regions not well covered by the drifters such as the southwest 

Gulf.

5.5 Verification with Oil Spill Trajectories

The intended use of this dataset is for oil spill risk analysis. Therefore, we ran simple 

experiments where we launched oil particles at the DwH oil well head at 1500 m and tracked 

them using the hindcast circulation fields for a month. The oil particles are characterized by a 

density and size, but are otherwise passive. About 10 particles were launched every 30 min for a 

duration of 1 month starting on April 20, 2010. Figure 14 shows the location of the particles after 

one month. The red line is the oil extent as analyzed by NOAA. We also show a satellite image 

of the oil obtained from MODIS (Moderate Resolution Imaging Spectroradiometer) Terra 

satellite on May 17. The image clearly indicates the model’s skill in capturing the overall shape 

and the southeast extension of the oil from the DwH location. The modeling of the surface oil 

slick is complicated and numerous forcing factors are involved transporting the oil at the surface. 

The particle distributions shown here are similar to the results presented in Mariano et al. (2011)

and LeHenaff et al. (2013) where the US Navy’s GOM-HYCOM velocity fields were used. 

6. Discussion and Conclusions

The aim of this work was to generate a 10-year (2003-2012) hindcast of the ocean state for the 

GoM for use in oil spill risk analysis. We have used the latest version of the HYbrid Coordinate 

Ocean Model (HYCOM)together with a multi-scale ensemble based assimilation scheme to 

constrain the model to observations by updating the model state once daily. The primary 

observation dataset used was along-track sea level anomalies. This was complemented by 

gridded sea surface temperature, profiles of salinity and temperature and climatology profiles. 

An assessment of the data assimilative hindcast shows robust model skill in reproducing the 

various features of the circulation in the GoM. There is a good agreement between the hindcast 

and independent estimates of Loop Current behavior inferred from altimetry and ocean color 

both in terms of temporal and spatial characteristics. Furthermore, the assimilation of a limited 

number of in-situ and climatological profiles of temperature and salinity is effective in 

controlling bias and reduces errors in the vertical density structure in the model. The upper ocean 

temperature and salinity match observation to within 1.5°C and 0.3 psu both well below the 

climatological variance. These statistics are calculated using all of the observations which are

much more in number than the actual number assimilated and thus constitute an robust check on 

the performance of the hindcast. Comparisons with independent drifter derived velocity data 

indicate significant correlation (>0.6) between the drifter data and model velocities. Simplified 

simulations of the DwH event captures the main features of surface plume indirectly verifying 

the underlying model circulation fields. We have not evaluated the deeper circulation as much as 

the upper ocean due to a lack of data. Velocity data from the deeper portion of the Gulf will be 

useful for such evaluations.
18

Figure 14. Location of the particles one month after release at the wellhead.

A satellite image of the oil obtained from MODIS Terra satellite on May 17 (top 

left). The hindcast is in the top right panel with the red line being the oil extent 

as analyzed by NOAA. The global HYCOM reanalysis and the GOM-HYCOM 

hindcast run by NRL are also shown for comparison. 

7. Project Deliverables

The following products have been prepared as part of the deliverable for this project and can be 

copied and delivered in a hard disk to BOEM. These deliverables will also be available at 

http://stargazer.coaps.fsu.edu/boem_gom_hindcast:
19

• The complete outputs of the hindcast dataset discussed including all variables in the 

format requested by the program officers.

• Model code and documentation along with a ready-to-run test case for 2010.
20

References

Barnier, B. Siefridt, L., Marchesiello. 1995. Thermal forcing for a global ocean model using a 

three year climatology of ECMWF analysis. J. Mar. Sys., 6, 363-380.

Bennett, A. F., 1992. Inverse Methods in Physical Oceanography. Cambridge, University Press, 

England.

Bleck, R., 2002. An oceanic general circulation model framed in hybrid isopycnal-cartesian 

coordinates. Ocean Modelling, 4, 55-88.

Carnes, R. M., 2009. Description and evaluation of GDEM-V 3.0. NRL Tech Report: 

NRL/MR/7330--09-9165.

Chassignet, E., P., Smith, L. T., Halliwell, G. R., Bleck, R., 2003. North Atlantic simulations 

with the hybrid coordinate ocean model(hycom): Impact of vertical coordinate choice, 

reference pressure and thermobaricity. J. Phys. Oceanogr., 33, 2504--22526.

Chassignet, E., P., Hurlburt, H., E., Smedstad, O., M., Barron, C., N., Ko, D., S., Rhodes, R., C., 

Shriver, J., F., Wallcraft, A., J., and Arnone, R., A., 2005. Assessment of Data assimilative 

ocean model in the Gulf of Mexico using ocean color. In Circulation in the Gulf of Mexico, 

Observations and Models, AGU Monograph Series, editors: Sturges and Lugo-Fernandez 

Chassignet, E., P., Hurlburt, H., E., Smedstad, O., M., Halliwell, G., R., Wallcraft, A., J., 

Metzger, E., J., B., O. Blanton, Lozano. C., Rao, D., B., Hogan, P., J., Srinivasan, A., 2006. 

Generalized vertical coordinates for eddy-resolving global and coastal ocean forecasts. 

Oceanography, 19, 22-33.

Cooper, M. and K. Haines, 1996. Altimetric assimilation with water property conservation. J.

Geophys. Res , 101, 1059-1078. 

Gaspari, G. and Cohn, S.E., 1999. Construction of correlation functions in two and three 

dimensions. Quart.J.Royal Met. Soc., 125, 723–757.

Halliwell, Jr., G. R., A. Srinivasan, V. Kourafalou, H Yang, D. Willey, M. Le Henaff, and R. 

Atlas, 2014. Rigorous evaluation of a fraternal twin ocean OSSE system in the open Gulf of 

Mexico. J. Atmos. Ocean. Technol. DOI: 10.1175/JTECH-D-13-00011.

Kara, A.B., Wallcraft, A.J., Martin, P. and Pauley R., 2009. Optimizing surface winds using 

QuikSCAT measurements in the Mediterranean Sea during 2000–2006. J. Mar. Sys., 78,

S119–S13.

Kantha, L., Choi, J., Schaudt, K., and Cooper, C.K., 2005. A regional data assimilative model for 

the Gulf of Mexico. In Circulation in the Gulf of Mexico, Observations and Models, AGU 

Monograph Series, editors: Sturges and Lugo-Fernandez

Leben, R.R., G.H. Born, and B.R. Engebreth , 2002. Operational altimeter data processing for 

mesoscale monitoring. Marine Geodesy, 25, 3-18.

Le Hénaff, M., Kourafalou, V.H., Paris, C., B., Helgers, J., Aman, Z., Hogan, P.J., and 

Srinivasan, A., 2013. Surface Evolution of the Deepwater Horizon Oil Spill Patch: Combined 

Effects of Circulation and Wind-Induced Drift. Environ. Sci. Technol., doi: 0.1021/es301570

Mariano, A.J., V.H. Kourafalou, A. Srinivasan, H. Kang, G.R. Halliwell, E.H. Ryan, M. Roffer 

2011. On the modeling of the 2010 Gulf of Mexico Oil spill. Dynamic of Ocean and 

Atmospheres, 52, 322–340

Mesinger, F., DiMego, G., Kalnay, E., Mitchell, K., Shafran, P. C., Ebisuzaki, W., and Shi, W. 

2006. North American regional reanalysis. Bulletin of the American Meteorological Society, 

87(3), 343-360.
21

Oey, L.-Y., Ezer, T., and Lee, H.C., 2005, Loop current Rings and Related Circulation in the 

Gulf of Mexico: A Review of Numerical models and future challenges. In Circulation in the 

Gulf of Mexico, Observations and Models, AGU Monograph Series, editors: Sturges and 

Lugo-Fernandez 

Srinivasan, A., E.P. Chassignet, L. Bertino, J.M. Brankart, P. Brasseur, T.M. Chin, F. Counillon, 

J.A. Cummings, A.J. Mariano, O.M. Smedstad, and W.C. Thacker, 2011. A comparison of 

sequential assimilation schemes for ocean prediction with the HYbrid Coordinate Ocean 

Model (HYCOM): Twin experiments with static error covariance. Ocean Modelling, 37, 

doi:10.1016/j.ocemod.2011.01.006, 85-111.

Velissariou, P., 2011. Development of a high resolution bathymetry for the Gulf of Mexico, 2011 

available from ftp://coaps.fsu.edu/pub/takis

Wunsch, C., 1996. The Ocean Circulation Inverse Problem. Cambridge University Press, 

England.
22

Appendix A. Optimizing the Atmospheric Forcing

Near-surface variables from atmospheric reanalysis – wind speed, air temperature, 

humidity, downward radiation and precipitation – are used to force ocean models. Using 

the output of atmospheric reanalysis allows the use of high-temporal-resolution winds and 

the specification of consistent fluxes since all of these are obtained through the surface 

boundary layer parameterization used in the atmospheric model. However, these variables 

are characterized by significant biases and errors that can lead to problems when used for 

forcing ocean models in long-term hindcasts or reanalysis. Regular assimilation of 

observations can compensate for these errors to some extent but the long time-scale errors 

associated with the biases in atmospheric forcing in general have different spatial 

correlations than the typical forecast errors that the assimilation schemes are designed to 

correct. This eventually leads to inconsistency between the assimilated solution and the 

model state resulting in a discontinuous model trajectory during long integrations. 

Although errors are controlled, they are typically greater than what would be expected 

using an unbiased model and when all of the information present in the observations is 

retained. In general, it is better to identify and correct the causes of bias in model before 

assimilating data as this leads to a smoother model trajectory during a long hindcast or 

reanalysis. The aim is to correct the source of errors rather than the consequences of the 

errors. By correcting the forcing functions, the problem of long-term trends and biases can 

be addressed in a manner akin to optimal control methods while ensuring the internal 

consistency of the corrections with the original reanalysis products.

A-1. Radiation Flux Corrections for Thermal Forcing

Excessive values of net solar radiation in the forcing data or incorrect vertical mixing 

parameterization can cause overheating of the model ocean. We assume that the 

overheating is mainly due to incorrect surface fluxes because the vertical parameterization 

used here has been used previously without the deleterious effects. Therefore, we estimate 

a constant in time, radiation flux correction using a method similar to the one described by 

Barnier et al. (1995). These corrections are obtained from the model free run climatological 

SSTs and NOAA’s optimum surface temperature analysis. The procedure involves several 

steps and is described below: 

• Bin average HYCOM SST for each month 1-12 for years 2006 and 2007. Size of the 

bin is 0.25 deg consistent with NARR and AVHRR resolution. The years 2006-2007 

are chosen to form the climatological errors after allowing for a 3 year model 

spinup.

• From monthly means of binned model SSTs 

• Compute the mean (over the 12 months) error of SST with respect to corresponding 

monthly mean AVHRR SST. 

• Convert the mean surface temperature error to surface flux correction is given by:

Q = ρ*Cp* Z * ΔSST⁄dt

In the above formula: 
23

ρ =1025 kg/m3

Cp~4000 J/(kg deg),

Z, is a typical thickness where the heat flux is distributed, 30 meters 

ΔSST is a monthly (mean value) of the error. 

Assuming that the flux acts during one month, time=2.6×106 sec, it would produce a 

temperature increase ΔT (the error) that we want to compensate. 

 Q = (1025 kg/m3 x 4000 (J/kg degK) x 30 m / 2.6×106s) ΔSST

Thus we get a factor of 47 Watts/m2 degK. The sign is negative because the error is SST 

observed – SST model. As can be seen the model is warmer than the AVHRR surface 

analysis during all times of the year and in most of domain. The climatological error is least 

about 1 deg during the winter times and peaks during the summer to about 1.8 deg. The 

mean error over year is approximately 1.25 degree, which translates into a net radiation 

flux correction of - 55.0-60W/m2. This correction is added to the total (downward) 

radiative surface flux: Oshortwv + Qlongwv. The magnitude of the correction estimated 

here is similar to the corrections estimated for ECMWF and NCEP reanalysis products. The 

correction term is constant in time and consistent with Barnier et al. (1995) who show that 

in areas not influenced by the monsoon regime, the pattern of the flux correction does not 

show drastic seasonal variations, and the pattern of the various monthly means globally 

follows that of the annual mean.

A-2. Wind Speed Correction Based on Linear Regression to 

QuickSCAT Winds

Winds just above the sea surface play an important role in many features of the Gulf of 

Mexico. Ocean general circulation models (OGCMs) generally use wind forcing to simulate 

upper- and deep-ocean fields, including sea-surface temperature (SST) and surface 

circulation. For applications such as oil spill risk analysis it is important to have accurate 

wind speed and direction at the sea surface (i.e., at 10-m height) to force the ocean model. 

Wind measurements with finer spatial resolution from the SeaWinds scatterometer sensor 

on the Quick Scatterometer (QuikSCAT) satellite over the global ocean are available since 

July 1999. However, the twice-daily QSCAT winds cannot be used directly for forcing the 

model since there can be data voids, depending on the coverage of the satellite passes. 

However, a correction based on the QSCAT monthly means can improve the accuracy of the 

3-hourly NARR winds, which can then be used for forcing. Therefore we optimize the hightemporal-resolution winds from the NARR reanalysis dataset using the (QuikSCAT) 

satellite data using a procedure detailed in Kara et al. (2009). The purpose of this 

optimization is to reduce possible errors in the NARR winds under the assumption that the 

QSCAT winds represent the truth. A comparison of the mean QSCAT and NARR wind speeds 

and standard deviations for 2003-2007 show that NARR winds are weaker than the QSCAT 

winds by about 2 m/s. To form the optimized winds we perform a linear regression 

analysis between the monthly mean NARR and QSCAT winds at each grid point. We derive 

a corrected NARR winds as NARR_corrected = aQSCAT + b, where a is the slope, and b is the 

intercept in ms-1. This is based on the fact that these winds are highly correlated with 
24

correlation coefficients of 0.85 based on a time series of 60 monthly mean winds almost 

everywhere over the Gulf of Mexico from 2003 through 2007. The slope and intercept 

values for the least squares lines are then computed at each grid point. We derive the 

corrected NARR wind speed at each grid point using the slope and intercept of the 

regression line. Figure 15 is a scatter plot between QSCAT and the NARR uncorrected and 

adjusted NARR winds. There are a total of 60 monthly mean wind values during the fiveyear time period. Significant improvement is evident from the adjusted NARR winds, 

resulting in values almost identical to the QSCAT winds. In fact, the mean bias between the 

basin-averaged QSCAT and NARR winds is zero during 2003–2007.

 

Figure 15. Scatter plot between monthly mean QSCAT winds and the NARR 

winds.

The adjusted NARR winds show good agreement with QSCAT winds with 

negligible bias.

We have chosen to compute the optimizations to the surface forcing based on 2003-2007 free run 

of the model. Although this time frame does not cover the entire hindcast time frame, we believe 

that it is sufficient to identify any existing systematic biases. We expect that any remaining errors 

will be further corrected by data assimilation.
25

Appendix B. Observation Error Specification and processing

Data assimilation typically involves the calculation of an analyzed state vector xa, by 

combining a background field x and a vector of observations d. The optimal combination is 

obtained by weighting the background field and observations with the inverse of their 

error variances, respectively. Therefore every observation that is assimilated must be 

assigned an error variance. Observation errors generally include measurement errors, 

representation errors and an age error. The observation error variance used in the 

assimilation step is then specified as the sum of the variances of various error components:

σooo

2 = σii iiiie 2 + σ r r r r

2 + σ a

2

Instrument measurement errors are dependent on the precision and biases of the 

instrument platforms. They are generally obtained by calibration of the instruments and 

are provided by the manufacturer. The instrument error standard deviations for different 

types of observations are listed in the table below. 

Table 4. Summary of the assumed standard deviations of instruments errors.

Instrument Instrument Error Standard Deviation

CTD, Float, Glider temperature 0.05° C

CTD, Float, Glider, salinity 0.02 PSU

XBT Temperature 0.2° C

SLA –TP/J1/J2 0.03 m

Envisat, Cryosat 0.05 m

Representation errors arise due to the fact that measurements contain information on all 

space and time scales while the model only represents a finite number of scales. Therefore 

these errors are representative of sub-grid scale phenomena and are a function of model 

resolution. Their magnitude will generally decrease as model resolution increases. Age 

errors arise in cases where measurements made within a time window are used for the 

analysis. In this situation many of the observations will be collected at times different from 

the assimilation time. Observations collected at times closer to the analysis time should 

have greater weight in the analysis procedure. This is incorporated into the assimilation 

methodology by using the absolute value of the difference in time between the 

measurement and the analysis time as age error. The variances of the observations are then 

increased in some specified manner as the age of the observations increase. 

As a first step in the observation preprocessing cycle, we check and flag any obvious 

outliers with values significantly different from the model and observation climatology. 

Data that lies outside five standard deviations from the intra-seasonal model anomalies are 

excluded. After this step, the instrument, representation and age errors for each type of 

observations are estimated separately and then combined for the assimilation step. For our 
26

high-resolution model (3.5 km) all observing platforms have a footprint that is coarser than 

the model resolution. We assume that the observations represent some average of the 

region of the swath or track and average the model to match the footprint of the 

observations. 

The in-situ data sampling is highly irregular. Many of the profiles were obtained were extremely 

close in both space (less than a 0.1°) and time (less than 1day). Further almost half of the profiles 

only extend to the upper 50 meters of the water column. Although all profiles are flagged with 

data quality markers but many of the data quality flags were erroneous and a large amount of 

time was spent in deciding which profiles to assimilate and which one to leave out. In light of the 

above issues, a conservative thinning procedure used in the pre-processing the data. We 

assimilate no more than one profile for a 0.25° grid box. When more than one profile per grid 

box is present, we assimilated the closest and the one with most data points (with better vertical 

sampling). Because our correlations are derived from the model and the ensemble size is limited 

there is a doubt on the robustness of the correlations in the shallower regions of the domain. So 

we decided to limit the profile assimilation to areas where water depth is at least 50 meters or 

greater. 
27

Appendix C: Tendral-Statistical Interpolation system (T-SIS)

implementation details for HYCOM

This appendix details the implementation details of the Tendral-Statistical Interpolation 

system for HYCOM. T-SIS is a statistical interpolation package for use with ocean circulation 

models in analysis, forecasting and system evaluation applications. T-SIS can be used with all 

ocean model types, and the version used herein is optimized for assimilation into the hybrid 

Lagrangian vertical coordinate layers of HYCOM. Figure 16 shows the schematic of the 

statistical interpolation algorithm used in the T-SIS assimilation scheme.

In general, the model state vector used in estimation procedures contains all of the prognostic 

variables. However, the state vector used for T-SIS is a subset of the HYCOM prognostic 

variables, specifically layer thickness, layer temperature, layer salinity, layer density, and the 

diagnosed sea-surface height (SSH) anomaly. In addition, the state vector is further subdivided 

into three sub-vectors, one consisting of SSH anomaly, another of layer thickness, and another 

consisting of layer temperature, salinity, and density. Each sub-vector is assumed to be 

uncorrelated with the others, making the covariance matrix block diagonal. 

The SSH anomaly field is not directly assimilated because in HYCOM, it is diagnosed from the 

prognostic bottom pressure and internal density fields. Instead, a layerized version of the Cooper 

and Haines (1996) procedure is used to adjust model layer thicknesses in the isopycniccoordinate interior in response to SSH anomaly innovations. Prior to calculating SSH 

innovations, the mean dynamic topography is added back into the altimetry observations. 

To optimize system performance for the HYCOM Lagrangian vertical coordinate system 

(essentially a stack of shallow water layers), subsurface profile observations are first layerized 

(re-mapped onto the model hybrid isopycnic-sigma-z vertical coordinate system) prior to 

assimilation. The analysis procedure then updates each layer separately in a vertically decoupled 

manner. For temperature profiles that do not have corresponding salinity profiles, synthetic 

salinity profiles are generated from climatological T-S relationships to permit layerization. The 

situation with velocity components is more complicated given that they are decomposed into 

barotropic and baroclinic components in the model. The barotropic velocity components are not 

included in the estimation procedure because the observation types, sampling frequencies, and 

assimilation time window are not appropriate to constrain barotropic velocity. In this initial 

version of T-SIS, baroclinic velocity components are also excluded from the estimation 

procedure since the cross-correlations required to update them are typically not robust enough. 

Instead, a geostrophic velocity update increment is calculated from layer pressure increments as 

a post-processing step after all other fields have been updated, and is then used to adjust both the 

barotropic and baroclinic velocity components.

The above modifications do make the estimation procedure less than fully multivariate, but are 

effective in the absence of robust cross-correlations between the sub-vectors. The analysis is 

performed in the model grid space using a simple, non-adaptive, distance-based localization 

using observations within the localization window around a particular grid point. A quasiGaussian, isotropic, distance-dependent localization function (Gaspari and Cohn, 1999) is used 

to impose a smooth localization of the error covariance and the innovations to yield a spatially 

continuous analysis. The localization radius, beyond which the ensemble-based covariance 
28

between two points is artificially reduced to zero, is uniform in space and set to 150 km. This 

corresponds to an e-folding radius of about 90 km. If instances of negative layer thickness occur 

after performing an analysis cycle, they are corrected as a post-processing step.

Figure 16. Schematic of the statistical interpolation algorithm used in the TSIS assimilation scheme.

T-SIS offers the flexibility to be fully multivariate or use sub vectors for corrections. In the 

absence of robust cross correlations the full state vector is split into sub vectors. In the model’s 

pressure layers, T&S are corrected. In isopycnal layers, thickness and one of either T/S is 

corrected. Profile data is first assimilated to obtain an intermediate state is then corrected with 

altimeter data. Velocity is not directly updated instead geostrophic velocity corrections are 

estimated from corrected state. The analysis process is designed to ensure that the mean 

model state is constrained mainly by T/S profiles while the variability around the mean is 

then constrained by the altimeter data.
29

Appendix D: System Documentation

D-1. Overview

The Gulf of Mexico hindcast system configured as a ready-to-run test case for the year 2010 is provided 

as a project deliverable. The test case for 2010 will be shipped to BOEM in a hard disk along with these 

instructions and the complete 2003-2012 hindcast dataset. The website: 

http://stargazer.coaps.fsu.edu/boem_gom_hindcast contains a copy of the package shipped to BOEM. All 

updates will be available through this site.

D-2. Contents of the Package

The package as distributed contains the hindcast datsaset, model, data assimilation codes 

observations and other input files needed to runs a 2010 Gulf of Mexico hindcast and animations 

and data used for validations. The contents of the package are distributed in several folders as 

below:

bin -- executables produced after compiling the sources

clim -- climatology for model runs and assimilation

data -- hindcast data ordered by year 

force -- prebuilt atmospheric forcing 

nlists -- input files for running the system

nest -- lateral boundary forcing files for 2010

obs -- observations for assimilation

scripts -- python scripts to run the hindcast

src -- HYCOM and TSIS source codes

restart -- model initial conditions for Jan 1, 2010

topo -- contains all model grid definition files

validation – animations and other data used for validation

D-3. Compiling the package

D-3.1 Required Compilers and Libraries

The package requires a modern fortran compiler and support libraries to compile the sources. 

Before compiling the code on your computer, check to see if you have a working modern fortran 

(90 or later compiler), netCDF and Message Passing Interface (MPI) libraries installed (for 

parallel processing). The code has been tested with Intel, Portland Group and gfortran compilers. 

It might require tweaks for other compilers. If netCDF and MPI are installed in a non-standard 

location (e.g. not in /usr/local/), then the paths to netCDF library, and to its include/ directory 

must be provided. 
30

D-3.2 Compiling the HYCOM ocean model

The first step in compiling the package is to compile the ocean model. The model sources are in 

the src/model folder. HYCOM requires that the number of processors and domain partitioning 

information be specified at compile time in the dimension.h file. By default the package comes 

configured for running on 32 processors using Intel compilers. The details of the compilers and 

libraries are included in the Makefile headers. For Intel compilers the file intelIFC in the source 

folder contains the compiler and library information as shown below:

#

# ---------------------------------------------------------------------

# common definitions for Intel Linux/IFC, MPI, real*8

# --------------------------------------------------------------------

# MACROS DESCRIPTIONS:

#

# FC: Fortran 90 compiler.

# FCFFLAGS: Fortran 90 compilation flags.

# CC: C compiler.

# CCFLAGS: C compilation flags.

# CPP: cpp preprocessor (may be implied by FC).

# CPPFLAGS: cpp -D macro flags.

# LD: Loader.

# LDFLAGS: Loader flags.

# EXTRALIBS: Extra local libraries (if any).

#

FC = mpif90

FCFFLAGS = -vec_report0 -O3 -msse3 -r8 -i_dynamic -w

CC = gcc

CCFLAGS = -O

CPP = cpp -P

CPPFLAGS = -DREAL8 -DMPI -DSERIAL_IO -DTIMER -DENDIAN_IO -I/opt/netcdf/include

LD = $(FC)

LDFLAGS = $(FCFFLAGS)

EXTRALIBS =-L/opt/netcdf/lib -lnetcdff -lnetcdf

Edit the above file to reflect the local environment. Once this is done the model can be compiled 

as is by the following commands as shown in the example below:

ashwanth@n1:/raid/FS/GOMFS> cd src/model

ashwanth@n1:/raid/FS/GOMFS/src/model> make

This should produce the executable hycom. Once the executable is compiled it should be 

moved to the bin directory.

ashwanth@n1:/raid/FS/GOMFS/src/model> mv hycom ../../bin/.
31

D-3.3 Compiling the T-SIS Assimilation Codes

The next step in compiling the package is to compile the T-SIS codes. The sources are in the 

src/tsis folder. T-SIS is designed to be compiled once to run on any number of processors. By 

default the package comes configured for compiling using Intel compilers. The details are in the 

Makefile headers in the arch directory. The T-SIS code requires LAPACK (www.netlib.org). 

The code comes configured for using the Intel Math kernel Libraries. These have to be edited to 

suit the local environment

ashwanth@n1:/raid/FS/GOMFS/> cd src/tsis/arch

For Intel compilers the file ifort_x86_64_mpi in the source folder contains the compiler and 

library information as shown below:

F90 = mpif90

F77 = mpif90

CC = gcc

CPP = cpp -P

LD = $(F90)

INC = -I/usr/local/include -I/opt/intel/mkl/include 

LIBS = -L/usr/local/netcdf/lib64 -L/opt/intel/mkl/lib -lnetcdf -lnetcdff -lmkl_intel_lp64 -

lmkl_sequential -lmkl_core 

FCFLAGS = -g -O2 -convert big_endian -w $(INC)

CCFLAGS = -O2

CPPFLAGS = -DMPI

LDFLAGS = $(FCFLAGS)

Once the arch information is edited to reflect the local environment, the source is ready to be 

compiled as below:

ashwanth@n1:/raid/FS/GOMFS/src/tsis/src/> make

Once the code is compiled the executables “tsisx”, “setobs and “setclim” should be moved to the 

bin directory

D-4. Setting up the 2010 Hindcast

The scripts directory contains several python scripts that glue the various components and moves 

data between the different components to produce a seamlessly functioning system. The only 

files that a user has to edit are the defs.py and the hindcast.py files. The defs.py file contains 

information on the top level the directory structure as shown below:

# defs.py defining filesystem locations

# version: 0.1

# July 2012

# Report bugs to: a.srinivasan@tendral.com

#############################################################

#Modify locations below to match your settings #

#############################################################
32

TLDIR="/raid/FS/GOMFS"

BIN="/raid/FS/GOMFS/bin"

dbnlst="/raid/FS/GOMFS/nlists"

dbclim="/raid/FS/GOMFS/clim"

dbnst="/raid/FS/GOMFS/nest"

dbfrc="/raid/FS/GOMFS/force"

mgfile="/raid/FS/GOMFS/topo/gridinfo.nc"

The above file has to be edited to specify the top-level directory. For example if the package is 

unpacked in /home/john, TLDIR should be specified as: TLDIR="/home/john/GOMFS. Similarly 

the other locations should be specified to match the local settings. 

The next file to edit is hindcast.py which controls the entire system:

#!/usr/bin/python

# hindcast.py driver script to control the assimilation system

# version: 0.1

# July 2012

# Report bugs to: a.srinivasan@tendral.com

import sys

import datetime as dt

from defs import *

from limits import dolimits

from hycom import dohycom

from clean import do_hindcast_clean

from assim import doprep,doanal

from arch2nc import doarch2nc

format = '%Y%m%d'

#############################################################

#Modify start date no of days and the number of processors #

#############################################################

dtgs=dt.datetime.strptime(str(20100101), format)

ndays=365

npes_hycom=32

npes_tsis=32

#############################################################

# End of user modifications #

#############################################################

dtgref=dt.datetime.strptime(str(19001231), format)

jdateDiff = dtgs - dtgref

je=jdateDiff.days

jd_start=je

adtg=dtgs
33

for day in range(jd_start, jd_start+ndays): 

 start=day

 end=start+1

 dolimits(scr_dir,start,end)

 dohycom(scr_dir,npes_hycom)

# do assimilation 

 adtg=adtg+dt.timedelta(days=1)

 sadtg=adtg.strftime('%Y%m%d')

 doprep(scr_dir,sadtg)

 doanal(scr_dir,npes_rls)

 do_hindcast_clean(scr_dir)

# convert outputs to diagnostic netcdf format

doarch2nc(scr_dir,cdtg)

sys.exit(0)

The only lines requiring user modifications are the start date and the number of processors 

listed below:

dtgs=dt.datetime.strptime(str(20100101), format)

ndays=365

npes_hycom=32

npes_tsis=32

Once this is set the system is ready for running and can be launched as

ashwanth@n1:/raid/FS/GOMFS/scripts> ./hindcast.py >& hindcast.log &

That’s it! This should run the hindcast for 2010. A scratch folder will be created where the model 

will run and produce outputs. These can be moved to a location of your choice.

D-5. Reference Model Inputs Parameters

The HYCOM model uses a text file called “blkdat.input” which contains the input parameters. 

Normally the default file in the package should be sufficient but if modifications are necessary 

then the HYCOM documentation at (www.hycom.org) should be consulted. The blkdat.input file 

used in the hindcast is reproduced below.

D-5.1 HYCOM Reference blkdat.input

#NARR forcing; LWcorr; precip; SSSrlx;

#Sigma0; Nested in ATLd0.08 7.2; KPP mixed layer; SeaWiFS mon KPAR;

#S-Z(15-5): dp00/f/x/i=3m/1.125/12m/1m; ds=1m/1.125/4m; src_2.2.03;

#12345678901234567890123456789012345678901234567890123456789012345678901234567890
34

# 22 'iversn' = hycom version number x10

# 008 'iexpt ' = experiment number x10

# 541 'idm ' = longitudinal array size

# 385 'jdm ' = latitudinal array size

# 400 'itest ' = grid point where detailed diagnostics are desired

# 200 'jtest ' = grid point where detailed diagnostics are desired

# 20 'kdm ' = number of layers

# 20 'nhybrd' = number of hybrid levels (0=all isopycnal)

# 15 'nsigma' = number of sigma levels (nhybrd-nsigma z-levels)

# 3.0 'dp00 ' = deep z-level spacing minimum thickness (m)

# 12.0 'dp00x ' = deep z-level spacing maximum thickness (m)

# 1.125 'dp00f ' = deep z-level spacing stretching factor 

# 1.0 'ds00 ' = shallow z-level spacing minimum thickness (m)

# 4.0 'ds00x ' = shallow z-level spacing maximum thickness (m)

# 1.125 'ds00f ' = shallow z-level spacing stretching factor 

# 1.0 'dp00i ' = deep iso-pycnal spacing minimum thickness (m)

# 6.0 'isotop' = shallowest depth for isopycnal layers (m), 

# 35.0 'saln0 ' = initial salinity value (psu

# 0 'locsig' = locally-referenced pot. density for stability (0=F,1=T)

# 0 'kapref' = thermobaric ref. state (-1=input,0=none,1,2,3=constant)

# 0 'thflag' = reference pressure flag (0=Sigma-0, 2=Sigma-2, 4=Sigma-4)

# 25.0 'thbase' = reference density (sigma units)

# 0 'vsigma' = spacially varying isopycnal target densities (0=F,1=T)

# 19.50 'sigma ' = layer 1 isopycnal target density (sigma units)

# 20.25 'sigma ' = layer 2 isopycnal target density (sigma units)

# 21.00 'sigma ' = layer 3 isopycnal target density (sigma units)

# 21.75 'sigma ' = layer 4 isopycnal target density (sigma units)

# 22.50 'sigma ' = layer 5 isopycnal target density (sigma units)

# 23.25 'sigma ' = layer 6 isopycnal target density (sigma units)

# 24.00 'sigma ' = layer 7 isopycnal target density (sigma units)

# 24.70 'sigma ' = layer 8 isopycnal target density (sigma units)

# 25.28 'sigma ' = layer 9 isopycnal target density (sigma units)

# 25.77 'sigma ' = layer 10 isopycnal target density (sigma units)

# 26.18 'sigma ' = layer 11 isopycnal target density (sigma units)

# 26.52 'sigma ' = layer 12 isopycnal target density (sigma units)

# 26.80 'sigma ' = layer 13 isopycnal target density (sigma units)

# 27.03 'sigma ' = layer 14 isopycnal target density (sigma units)

# 27.22 'sigma ' = layer 15 isopycnal target density (sigma units)

# 27.38 'sigma ' = layer 16 isopycnal target density (sigma units)

# 27.52 'sigma ' = layer 17 isopycnal target density (sigma units)

# 27.64 'sigma ' = layer 18 isopycnal target density (sigma units)

# 27.74 'sigma ' = layer 19 isopycnal target density (sigma units)

# 27.82 'sigma ' = layer 20 isopycnal target density (sigma units)

# 3 'iniflg' = initial state flag (0=levl, 1=zonl, 2=clim, 3=archv)

# 1 'jerlv0' = initial jerlov water type (1 to 5; 0 to use KPAR)

# 3 'yrflag' = days in year flag (0=360, 1=366, 2=366J1, 3=actual)

# 0 'sshflg' = diagnostic SSH flag (0=SSH,1=SSH&stericSSH)

# 0.0 'dsurfq' = number of days between model diagnostics at the surface 0

# 0.125 'diagfq' = number of days between model diagnostics

# 0.0 'tilefq' = number of days between model diagnostics on selected tiles

# 0.0 'meanfq' = number of days between model diagnostics (time averaged)

# 30.0 'rstrfq' = number of days between model restart output

# 15.21875 'bnstfq' = number of days between baro nesting archive input

# 15.21875 'nestfq' = number of days between 3-d nesting archive input

# 0.0 'cplifq' = number of days (or time steps) between sea ice coupling

#180.0 'baclin' = baroclinic time step (seconds), int. divisor of 86400

# 6.0 'batrop' = barotropic time step (seconds), int. div. of baclin/2

# 1 'incflg' = incremental update flag (0=no, 1=yes, 2=full-velocity)

# 90 'incstp' = no. timesteps for full update (1=direct insertion)

# 1 'incupf' = number of days of incremental updating input

# 0.125 'wbaro ' = barotropic time smoothing weight
35

# 0 'btrlfr' = leapfrog barotropic time step (0=F,1=T)

# 0 'btrmas' = barotropic is mass conserving (0=F,1=T)

# 16.0 'hybrlx' = HYBGEN: inverse relaxation coefficient (time steps)

# 0.01 'hybiso' = HYBGEN: Use PCM if layer is within hybiso of target density

# 3 'hybmap' = hybrid remapper flag (0=PCM, 1=PLM, 2=PPM, 3=WENO-like)

# 0 'hybflg' = hybrid generator flag (0=T&S, 1=th&S, 2=th&T)

# 0 'advflg' = thermal advection flag (0=T&S, 1=th&S, 2=th&T)

# 2 'advtyp' = scalar advection type (0=PCM,1=MPDATA,2=FCT2,4=FCT4)

# 2 'momtyp' = momentum advection type (2=2nd order, 4=4th order)

# -1.0 'slip ' = +1 for free-slip, -1 for non-slip boundary conditions

# 0.05 'visco2' = deformation-dependent Laplacian viscosity factor

# 0.0 'visco4' = deformation-dependent biharmonic viscosity factor

# 0.0 'facdf4' = speed-dependent biharmonic viscosity factor

# 0.0025 'veldf2' = diffusion velocity (m/s) for Laplacian momentum dissip.

# 0.01 'veldf4' = diffusion velocity (m/s) for biharmonic momentum dissip.

# 0.0 'thkdf2' = diffusion velocity (m/s) for Laplacian thickness diffus.

# 0.005 'thkdf4' = diffusion velocity (m/s) for biharmonic thickness diffus.

# 0.005 'temdf2' = diffusion velocity (m/s) for Laplacian temp/saln diffus.

# 1.0 'temdfc' = temp diffusion conservation (0.0,1.0 all dens,temp resp.)

# 2.e-5 'vertmx' = diffusion velocity (m/s) for momentum at MICOM M.L.base

# 0.05 'cbar ' = rms flow speed (m/s) for linear bottom friction

# 2.2e-3 'cb ' = coefficient of quadratic bottom friction

# 0.0 'drglim' = limiter for explicit friction (1.0 none, 0.0 implicit)

# 0.0 'drgscl' = scale factor for tidal drag (0.0 for no tidal drag)

# 500.0 'thkdrg' = thickness of bottom boundary layer for tidal drag (m)

# 10.0 'thkbot' = thickness of bottom boundary layer (m)

# 0.02 'sigjmp' = minimum density jump across interfaces (kg/m**3)

# 0.3 'tmljmp' = equivalent temperature jump across mixed-layer (degC)

# 30.0 'thkmls' = reference mixed-layer thickness for SSS relaxation (m)

# 0.0 'thkmlt' = reference mixed-layer thickness for SST relaxation (m)

# 6.0 'thkriv' = nominal thickness of river inflow (m)

# 20.0 'thkfrz' = maximum thickness of near-surface freezing zone (m)

# 0 'iceflg' = ice model flag (0=none,1=energy loan,2=coupled/esmf) 

# 0.0 'tfrz_0' = ENLN: ice melting point (degC) at S=0psu

# -0.054 'tfrz_s' = ENLN: gradient of ice melting point (degC/psu)

# 0.0 'ticegr' = ENLN: temp. grad. inside ice (deg/m); =0 use surtmp

# 0.5 'hicemn' = ENLN: minimum ice thickness (m)

# 10.0 'hicemx' = ENLN: maximum ice thickness (m)

# 0 'ntracr' = number of tracers (0=none,negative to initialize)

# 0 'trcflg' = tracer flags (one digit per tr, most sig. replicated)

# 480 'tsofrq' = number of time steps between anti-drift offset calcs

# 0.0 'tofset' = temperature anti-drift offset (degC/century)

# 0.0 'sofset' = salnity anti-drift offset (psu/century)

# 1 'mlflag' = mixed layer flag (0=none,1=KPP,2-3=KT,4=PWP,5=MY,6=GISS)

# 1 'pensol' = KT: activate penetrating solar rad. (0=F,1=T)

# 999.0 'dtrate' = KT: maximum permitted m.l. detrainment rate (m/day)

# 19.2 'thkmin' = KT/PWP: minimum mixed-layer thickness (m)

# 1 'dypflg' = KT/PWP: diapycnal mixing flag (0=none, 1=KPP, 2=explicit)

# 480 'mixfrq' = KT/PWP: number of time steps between diapycnal mix calcs

# 1.e-7 'diapyc' = KT/PWP: diapycnal diffusivity x buoyancy freq. (m**2/s**2)

# 0.25 'rigr ' = PWP: critical gradient richardson number

# 0.65 'ribc ' = PWP: critical bulk richardson number

# 0.7 'rinfty' = KPP: maximum gradient richardson number (shear inst.)

# 0.45 'ricr ' = KPP: critical bulk richardson number

# 0.0 'bldmin' = KPP: minimum surface boundary layer thickness (m)

#1200.0 'bldmax' = KPP: maximum surface boundary layer thickness (m)

# 0.7 'cekman' = KPP/KT: scale factor for Ekman depth

# 1.0 'cmonob' = KPP: scale factor for Monin-Obukov depth

# 0 'bblkpp' = KPP: activate bottom boundary layer (0=F,1=T)

# 1 'shinst' = KPP: activate shear instability mixing (0=F,1=T)
36

# 1 'dbdiff' = KPP: activate double diffusion mixing (0=F,1=T)

# 1 'nonloc' = KPP: activate nonlocal b. layer mixing (0=F,1=T)

# 0 'latdiw' = GISS: activate bot.enhan.int.wav mixing (0=F,1=T)

# 0 'botdiw' = GISS: activate bot.enhan.int.wav mixing (0=F,1=T)

# 0 'difout' = K-PROF: output visc/diff coffs in archive (0=F,1=T)

# 1 'difsmo' = K-PROF: activate horiz smooth diff coeffs (0=F,1=T)

# 50.e-4 'difm0 ' = KPP: max viscosity due to shear instability (m**2/s)

# 50.e-4 'difs0 ' = KPP: max diffusivity due to shear instability (m**2/s)

# 1.e-4 'difmiw' = KPP: background/internal wave viscosity (m**2/s)

# 1.e-5 'difsiw' = KPP: background/internal wave diffusivity (m**2/s)

# 10.e-4 'dsfmax' = KPP: salt fingering diffusivity factor (m**2/s)

# 1.9 'rrho0 ' = KPP: salt fingering rp=(alpha*delT)/(beta*delS)

# 98.96 'cs ' = KPP: value for nonlocal flux term

# 10.0 'cstar ' = KPP: value for nonlocal flux term

# 0.0 'cv ' = KPP: buoyancy frequency ratio (0.0 to use a fn. of N)

# 5.0 'c11 ' = KPP: value for turb velocity scale

# 2 'hblflg' = KPP: b. layer interp. flag (0=const.,1=linear,2=quad.)

# 2 'niter ' = KPP: iterations for semi-implicit soln. (2 recomended)

# 0 'fltflg' = FLOATS: synthetic float flag (0=no; 1=yes)

# 4 'nfladv' = FLOATS:advect every nfladv bacl. time steps (even, >=4)

# 1 'nflsam' = FLOATS: output (0=every nfladv steps; >0=no. of days)

# 0 'intpfl' = FLOATS: horiz. interp. (0=2nd order+n.n.; 1=n.n. only)

# 0 'iturbv' = FLOATS:add horiz. turb. advection velocity (0=no; 1=yes)

# 1 'ismpfl' = FLOATS: sample water properties at float (0=no; 1=yes)

#4.63e-6 'tbvar ' = FLOATS: horizontal turb. vel. variance scale (m**2/s**2)

# 0.4 'tdecri' = FLOATS: inverse decorrelation time scale (1/day)

# 2 'lbflag' = lateral barotropic bndy flag (0=none,1=port,2=input)

# 0 'tidflg' = TIDES: tidal forcing flag (0=none,1=open-bdy,2=bdy&body)

#00000000 'tidcon' = TIDES: 1 digit per (Q1K2P1N2O1K1S2M2), 0=off,1=on

# 0.06 'tidsal' = TIDES: scalar self attraction and loading factor

# 1 'tidgen' = TIDES: generic time (0=F,1=T)

#39813.00 'tidrmp' = TIDES: ramp time (days)

# 0. 'tid_t0' = TIDES: origin for ramp time (model day)

# 12 'clmflg' = climatology frequency flag (6=bimonthly, 12=monthly)

# 2 'wndflg' = wind stress input flag (0=none,1=u/v-grid,2,3=p-grid)

# 4 'ustflg' = ustar forcingflag (3=input,1,2=wndspd,4=stress)

# 4 'flxflg' = thermal forcing flag (0=none,3=net-flux,1,2,4=sst-based)

# 4 'empflg' = E-P forcing flag (0=none,3=net_E-P, 1,2,4=sst-bas_E)

# 0 'dswflg' = diurnal shortwave flag (0=none,1=daily to diurnal corr.)

# 1 'sssflg' = SSS relaxation flag (0=none,1=clim)

# 0 'lwflag' = longwave (SST) flag (0=none,1=clim,2=atmos)

# 0 'sstflg' = SST relaxation flag (0=none,1=clim,2=atmos,3=observed)

# 0 'icmflg' = ice mask flag (0=none,1=clim,2=atmos,3=observed)

# 0 'flxoff' = net flux offset flag (0=F,1=T)

# 0 'flxsmo' = smooth surface fluxes (0=T, 1=F)

# 0 'relax ' = activate lateral boundary nudging (0=F,1=T)

# 0 'trcrlx' = activate lat. bound. tracer nudging (0=F,1=T)

# 1 'priver' = rivers as a precipitation bogas (0=F,1=T)

# 0 'epmass' = treat evap-precip as a mass exchange (0=F,1=T)

D-5.2 T-SIS Reference Namelist

The assimilation component is controlled by a tsis.nlist name list. An annotated version is 

reproduced below for reference.
37

=====================T-SIS MASTER NAMELIST=========================

The model and first guess file name 

&model_data_type

model='hycom' 

FileName='restart_out' 

/

The observation location on disk

&obs_database_location

sla_data_location='/GOMFS/obs/altim/'

sst_data_location='/GOMFS/obs/sst'

insitu_data_location='/GOMFS/obs/insitu'

clim_data_location='/GOMFS/clim/'

fclim_profile_location='/GOMFS/clim/gom_gdem_profile_locations'

/

The observation window (can be different from model domain)

&obs_window ! 

lnmn=-97.5 ! min longitude 

lnmx=-77.0 ! max longitude

ltmn=18.10 ! min latitude

ltmx=31.0 ! max latitude 

/

The ensemble location and number of members

&ensnml

ensloc='/GOMFS/BASIS'

enslist='/GOMFS/BASIS/ens.list'

ntotalmembers=108

/

Three analysis parameters for each analysis variable: 

1)acdepth the critical depth in meters to stop assimilation when the local depth is less than 

this value

2)arlocal is the localization radius in kilometers 

3)acovscl is a scaling parameter [0,1] to strengthen or decrease the covariance between 

variables

&anlp

#! sla analysis parameters

 acdepth(1)=150.0

 arlocal(1)=150.0

 acovscl(1)=0.7
38

#! sst analysis parameters

 acdepth(2)=2.0

 arlocal(2)=50.0

 acovscl(2)=0.7

#! sss analysis parameters

 acdepth(3)=100.0

 arlocal(3)=150.0

 acovscl(3)=0.7

#! tem analysis parameters

 acdepth(4)=100.0

 arlocal(4)=150.0

 acovscl(4)=0.7

#! sal analysis parameters

 acdepth(5)=100.0

 arlocal(5)=150.0

 acovscl(5)=0.7

#! thk analysis parameters

 acdepth(6)=100.0

 arlocal(6)=150.0

 acovscl(6)=0.7

#! den analysis parameters

 acdepth(7)=100.0

 arlocal(7)=150.0

 acovscl(7)=0.7

#! uvl analysis parameters

 acdepth(8)=100.0

 arlocal(8)=150.0

 acovscl(8)=0.7

#! vvl analysis parameters

 acdepth(9)=100.0

 arlocal(9)=150.0

 acovscl(9)=0.7

/
The Department of the Interior Mission

As the Nation's principal conservation agency, the Department of the 

Interior has responsibility for most of our nationally owned public 

lands and natural resources. This includes fostering sound use of our 

land and water resources; protecting our fish, wildlife, and biological 

diversity; preserving the environmental and cultural values of our 

national parks and historical places; and providing for the enjoyment 

of life through outdoor recreation. The Department assesses our 

energy and mineral resources and works to ensure that their 

development is in the best interests of all our people by encouraging 

stewardship and citizen participation in their care. The Department 

also has a major responsibility for American Indian reservation 

communities and for people who live in island territories under US 

administration.

The Bureau of Ocean Energy Management

As a bureau of the Department of the Interior, the Bureau of Ocean 

Energy (BOEM) primary responsibilities are to manage the mineral 

resources located on the Nation's Outer Continental Shelf (OCS) in an 

environmentally sound and safe manner.

The BOEM Environmental Studies Program

The mission of the Environmental Studies Program (ESP) is to provide 

the information needed to predict, assess, and manage impacts from 

offshore energy and marine mineral exploration, development, and 

production activities on human, marine, and coastal environments.