Numerical Simulations and Observations of Surface Wave Fields under an Extreme Tropical Cyclone

The performance of the wave model WAVEWATCH III under a very strong, category 5, tropical cyclone wind forcing is investigated with different drag coefficient parameterizations and ocean current inputs. The model results are compared with field observations of the surface wave spectra from an airborne scanning radar altimeter, National Data Buoy Center (NDBC) time series, and satellite altimeter measurements in Hurricane Ivan (2004). The results suggest that the model with the original drag coefficient parameterization tends to overestimate the significant wave height and the dominant wavelength and produces a wave spectrum with narrower directional spreading. When an improved drag parameterization is introduced and the wave– current interaction is included, the model yields an improved forecast of significant wave height, but underestimates the dominant wavelength. When the hurricane moves over a preexisting mesoscale ocean feature, such as the Loop Current in the Gulf of Mexico or a warmand cold-core ring, the current associated with the feature can accelerate or decelerate the wave propagation and significantly modulate the wave spectrum.


Introduction
Tropical cyclone-generated wave fields are of interest both scientifically for understanding wind-wave interaction physics and operationally for predicting potentially hazardous conditions for ship navigation and coastal regions. There have been considerable efforts made to understand the characteristics of tropical cyclonegenerated surface waves through both measurements and numerical modeling. Wright et al. (2001) and Walsh et al. (2002) studied the spatial variation of hurricane directional wave spectra for both open ocean and landfall cases using the National Aeronautics and Space Admin-istration (NASA) Airborne Scanning Radar Altimeter (SRA). These measurements have provided detailed wave characteristics at a specific place and time. Moon et al. (2003) conducted a detailed comparison between WAVEWATCH III (WW3) wave model simulations and observations of the spatial distribution of hurricane directional wave spectra obtained from NASA SRA in Hurricane Bonnie (1998), a category 2-3 tropical cyclone on the Saffir-Simpson hurricane intensity scale (SSHS). Excluding shallow areas near the shore, the model results yielded good agreement with observations of directional spectrum as well as significant wave height, dominant wavelength, and dominant wave direction (wavelength and direction at the peak frequency of the wave spectrum). Later studies of Chao et al. (2005), Tolman and Alves (2005), and Tolman et al. (2005) found that WW3 overestimates the significant wave height under very high wind conditions in strong hurricanes. Moon et al. (2008) suggested that one of the reasons for the overestimation of the significant wave height is overestimation of the drag coefficient used in WW3 at very high winds. Comparing with buoy wave measurements during Hurricane Katrina in 2005 (SSHS category 5 in the Gulf of Mexico), they found that WW3 simulations with a reduced drag coefficient yielded more accurate simulations. During Hurricane Ivan (SSHS category 4-5 in the Caribbean Sea and Gulf of Mexico) in 2004, three sets of detailed SRA wave spectra measurements were collected by NASA through a joint effort between the NASA Goddard Space Flight Center and National Oceanic and Atmospheric Administration (NOAA)/ Atlantic Oceanographic and Meteorological Laboratory/ Hurricane Research Division (HRD). These observations, together with satellite measurements and National Data Buoy Center (NDBC) buoy time series, are used in this study to evaluate the WW3 predictions in extreme tropical cyclones. In particular, we investigate if an improved drag coefficient parameterization and inclusion of the effect of wave-current interaction may improve the wave predictions using WW3.
The outline of this paper is as follows: the wave and ocean models are described in section 2. The available observations are presented in section 3. Section 4 provides a new method for wind field specification. The results are presented and discussed in section 5. A summary and conclusions are given in section 6.

Methodology
a. The wave model The ocean surface wave model, WAVEWATCH III (Tolman 1998), has been used operationally at NOAA/ National Centers for Environmental Prediction (NCEP) for more than a decade. The model was validated over a global-scale wave forecast and a regional wave forecast (Tolman 1998(Tolman , 2002Tolman et al. 2002;Wingeart et al. 2001). WAVEWATCH III explicitly accounts for wind input, wave-wave interaction, and dissipation due to whitecapping and wave-bottom interaction, and solves the spectral action density balance equation for directional wavenumber spectra. The wave spectrum of the model is discretized using 24 directions and 40 intrinsic (relative) frequencies extending from 0.0285 to 1.1726 Hz, with a logarithmic increment of f n11 5 1.1f n , where f n is the nth frequency. The intrinsic frequency is related to the wavenumber (magnitude) k through the dispersion relation. The model domain is set to 58-328N in the latitudinal direction and 958-488W in the longitudinal direction, with a grid increment of 1 /128 in both directions for all experiments. Moon et al. (2004a,b) have shown that the drag coefficient used in WW3 is significantly larger than estimates based on a coupled wave-wind model (CWW) that explicitly integrates the waveform drag. The CWW results are more consistent with recent field and laboratory observations of the drag coefficient (Powell et al. 2003;Donelan et al. 2004;Black et al. 2007). Moon et al. (2008) found that using the CWW model in WW3 yielded improved wave predictions during Hurricane Katrina (2005).
In this study, we have conducted three sets of experiments with the WW3 model. In experiment A, the original drag coefficient parameterization in WW3 is used to force the wave model. In experiment B, the original drag coefficient parameterization has been replaced by the CWW model in calculating the wind input term in WW3. In the CWW model, WW3 is used to estimate the wave spectra near the peak. The spectra in the high-frequency range (equilibrium range) beyond the model resolution are parameterized by the analytical model of Hara and Belcher (2002). The resulting spectrum is then incorporated into the wave boundary layer model of Hara and Belcher (2004) to explicitly calculate the wave-induced stress vector, the mean wind profile, and the drag coefficient. The CWW model treats the wind stress as a vector quantity to consider the influence of dominant waves that propagate at a large angle to the local wind. It thus makes possible the estimation of the wind stress for any given surface wave field, even for the complex seas encountered under tropical cyclones.
In experiment C, we use the same setup for the wave model as in experiment B, but also introduce ocean currents that are produced by the ocean model described below in response to hurricane forcing. Funakoshi et al. (2008) also used a similar approach to study storm tides in the St. Johns River under Hurricane Floyd (1999). There are two significant ways the ocean current (U c ) impacts the wave field in the WW3 model. First, through the wind input term in the calculation of the wind stress. When ocean current is present, the 10-m wind velocity input (U 10 ) is replaced by the relative wind velocity U 10 -U c . Second, the wave action equation, which is solved in WW3, accounts for the modulation by the ocean current such that where N 5 c/v is the wave action spectrum, C g is group velocity vector, k is the wavenumber vector, u is the wave direction, s is a coordinate in the wave direction, m is the coordinate perpendicular to s, F represents all forcing terms, and U c is the ocean current at the depth of L/(4p) (L is the mean wavelength), which modifies the apparent phase speed of the wave train [please refer to Fan et al. (2009a)  b. The ocean model In this study, the ocean currents are calculated using the Princeton Ocean Model (POM). POM is a threedimensional, primitive equation model with complete thermohaline dynamics, sigma vertical coordinate system, and a free surface (Blumberg and Mellor 1987). We employed the version of POM used operationally in the Geophysical Fluid Dynamics Laboratory (GFDL)/ University of Rhode Island (URI) hurricane prediction system (Bender et al. 2007).
A realistic initialization of the 3D density and velocity fields in the ocean model is critical for proper simulation of the ocean response to a hurricane (Ginis 2002). The initialization method implemented in the GFDL/URI system is described in detail in Falkovich et al. (2005) and Yablonsky and Ginis (2008) and includes a realistic representation of the Loop Current (LC) and warm-and cold-core eddies in the Gulf of Mexico. In this method, the initial condition is first generated using the Generalized Digital Environmental Model (GDEM) monthly ocean temperature and salinity climatology (Teague et al. 1990), which has ¼8 horizontal grid spacing and 33 vertical z levels at depths ranging from 0 to 5500 m. Then, a feature-based modeling procedure (Yablonsky and Ginis 2008) is conducted to incorporate sea surface height anomaly (SSHA) observations. Positive SSHA features (also referred to as warm features) are regions where warm upper-ocean layer (warm ocean layer is often defined as from the surface down to the depth of the 268C isotherm) (Shay et al. 2000) is deeper than climatology. On the other hand, negative SSHA features are regions where the warm upper-ocean layer is shallower than climatology. These positive and negative SSHA features are frequently associated with warmand cold-core mesoscale eddies, respectively. In the feature-based procedure, multiple points along the LC path are specified, allowing the LC shape to be adjusted to match the observed shape derived from satellite altimetry. Then, the warm-and cold-core eddies in the Gulf of Mexico are incorporated by assuming they are elliptical in shape, with major and minor axes defined based on the SSHA data. For this study, we utilize the altimetry data from the Colorado Center for Astrodynamics Research (CCAR) Real-Time Altimetry Project through their Web site (http://argo.colorado.edu/ ;realtime/welcome/). The CCAR altimetry map on 12 September 2004 (shown in Fig. 15a) is used to ini-tialize the position and structure of the LC and a warmcore ring in the Gulf of Mexico (Fig. 15b).

Wave observations during Hurricane Ivan
Hurricane Ivan in 2004 was a classical, long-lived Cape Verde hurricane that reached SSHS category 5 strength three times. The hurricane track from 1200 UTC 8 September to 1200 UTC 16 September is shown in Fig. 1. Three sets of detailed SRA wave spectra measurements were collected by NASA through a joint effort between the NASA Goddard Space Flight Center and NOAA/HRD. The flight tracks of the NOAA aircraft carrying the SRA are shown in Fig. 1. Two sets of measurements were collected from 1615 to 2010 UTC 9 September and from 1040 to 1540 UTC 12 September when Ivan was crossing the Caribbean Sea and at its maximum intensity of category 5. The third set of measurements was done from 2030 UTC 14 September to 0330 UTC 15 September when Ivan entered the Gulf of Mexico.
The SRA measurements covered the region within about 28 of the hurricane eye. The SRA scanned a radar beam across the aircraft ground track to measure the elevation at 64 points on the sea surface. Sea surface topographic maps were produced from groups of SRA cross-track scan lines. The topography was then interpolated to a north-and east-oriented 256 3 256 rectangular grid of 7-m spacing centered on the data. The elevations in the uniform grid were transformed by a two-dimensional fast Fourier transform (FFT) with wavenumber spectral resolution of 0.0035 rad m 21 . The encounter wave spectra were Doppler corrected and the 1808 ambiguous spectral lobes were deleted. Wright et al. (2001) and Walsh et al. (2002) describe the process in detail.
Two satellites, Envisat-1 and ERS-2 (in the same orbit as Envisat-1 and trailing it by about 28 min), approached within about 90 km of the eye of Hurricane Ivan at 0338 UTC and 0406 UTC 15 September (Fig. 1). Both satellites carried radar altimeters that documented wind speed and wave height along their ground track. Also, five NDBC buoys, located within 48 of the hurricane track ( Fig. 1), documented significant wave height time series through the passage of the hurricane. All these data will be used to evaluate our model results.

Hurricane wind field specification
The wind fields during Hurricane Ivan are obtained from the NOAA/HRD real-time wind analysis (HWIND) and interpolated into 0.5-h intervals to input into both the WW3 and POM models. HWIND is an integrated tropical SEPTEMBER 2009 F A N E T A L .
cyclone observing system in which wind measurements from a variety of observation platforms are used to develop an objective analysis of the distribution of wind speeds in a hurricane (Powell et al. 1998). The wind data in gridded form are available at the HRD Web site for all hurricanes in the Atlantic basin since 1994 (www.aoml. noaa.gov/hrd/data_sub/wind.html). The HRD winds with the spatial resolution of about 6 km 3 6 km, covering an area of about 88 3 88 in latitude-longitude around the hurricane's center, are provided at intervals of every 3 or 6 h. This frequency is not sufficient to force a numerical model and therefore the wind data need to be interpolated in space and time. If a hurricane rapidly intensifies or its size undergoes significant changes in a short period of time, direct time/ space interpolation may result in distortion of wind fields. Any error in the input wind field will result in an error in the computed wave field because wind waves are very sensitive to small variations in the wind input. To illustrate this sensitivity, let us consider fully developed wave conditions in which the significant wave height H s is roughly proportional to the square of the wind speed and wave energy is roughly proportional to the wind speed cubed. A 10% bias in the surface wind speed may cause ;20% error in H s and ;35% in wave energy.
Here, we introduce a new interpolation method (hereafter called ''normalized interpolation'') of the HRD wind fields in time/space with minimum distortions of the hurricane wind field. For simplicity, we illustrate below the normalized interpolation method along one radial direction of the hurricane, which can easily be applied to a 2D hurricane wind field. Consider two radial profiles of wind speed W1(R1) and W2(R2) at two different times t 1 and t 2 with their maximum wind speed located at R1 and R2 correspondingly (Fig. 2a). A simple averaging of the two profiles at time (t 1 1 t 2 )/2 would result in the dashed line in Fig. 2a, which is clearly not a good approximation of the hurricane radial wind profile. In our method we first normalize the radial distance from the hurricane center by the radius of the maximum wind speed, so that in the normalized coordinate, both W1 and W2 have their maximum wind speed at the normalized distance 1.0 (Fig. 2b). If we interpolate these two wind profiles at time (t 1 1 t 2 )/2 (dashed line in Fig. 2b), the wind profile is not distorted like the one in Fig. 2a. Since the radius of the maximum wind speed for the interpolated wind profile is simply (R1 1 R2)/2, we use this radius to obtain the desired interpolated wind speed profile at time (t 1 1 t 2 )/2 in the dimensional coordinate as illustrated by the dashed line in Fig. 2c. An example of the 2D interpolated wind field at 1800 UTC 9 September and the HRD winds at 1330 and 1930 UTC 9 September is shown in Fig. 3.

a. Wave parameters
We first compare the H s , dominant wavelength (DWL), and dominant wave propagation direction (DWD) in WW3 with the SRA measurements. The H s is the standard output of the wave model. To obtain the DWL and DWD, the model directional frequency spectrum is transformed into the same wavenumber space as the SRA measurements using Jacobian transformation, and the location of the wave spectrum peak, which corresponds to the DWL and DWD, is determined using a parabolic interpolation. The DWL, DWD, and H s are interpolated from the uniform model grid in WW3 to the nonuniformly spaced SRA measurement locations using the cubic interpolation method. For each SRA location, these parameters are then linearly interpolated in time to obtain the DWL, DWD, and H s at the measurement time.
Comparisons between the model results in experiments A, B, and C and the observations at all the SRA measurement locations along the flight tracks are shown in Fig. 4 for the measurements on 9 September, Fig. 5 for the measurements on 12 September, and Fig. 6 for the measurements on 14-15 September. The model significant wave height and the dominant wavelength in experiments A, B, and C are plotted in Fig. 7 against the corresponding SRA data for all SRA measurement locations for the periods of 9, 12, and 14-15 September. The root-mean-square errors, defined as rmse 5 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1/Nå(x Model À x Observation ) 2 p , between the SRA measurements and the model results are presented in Table 1.
The model DWDs in all three experiments are very close to each other and match very well with the FIG. 2. Diagrams of wind profile interpolations. (a) Solid lines are wind profiles W1 and W2 vs radial distance at times t 1 and t 2 with maximum wind speed at R 1 and R 2 ; dashed line is wind profile obtained at time (t 1 1 t 2 )/2 using direct time/space interpolation. (b) Solid lines are wind profiles W1 and W2 vs radial distance normalized by the radius of maximum wind at times t 1 and t 2 ; dashed line is wind profile obtained at time (t 1 1 t 2 )/2 using normalized interpolation. (c) Solid lines are wind profiles W1 and W2 vs radial distance at times t 1 and t 2 with maximum wind speed at R 1 and R 2 ; dashed line is wind profile obtained at time (t 1 1 t 2 )/2 using normalized interpolation.
SRA observations during all three flight periods (Figs. 4b, 5b, and 6d). This indicates that the wind stress parameterization-based CWW model and wave-current interaction have negligible effects on the DWD predictions. The locations where the wave propagation directions differ by more than 108 are generally in the rear right quadrant of the hurricane where there can be two or three comparable peaks in the observed spectrum while the model spectrum has a smooth, one-peak structure.
From the H s comparison along the flight track on 9 September, we can see that at the locations in the rear quadrants of the hurricane, where the H s values are small (less than 5 m), the model predictions in experiment A agree very well with the SRA observations (Figs. 4d,7a). This is consistent with the fact that WW3 has been extensively calibrated and validated under low and moderate wind conditions. The H s values from experiments B and C are also very similar to those from experiment A at these locations (Fig. 4d), indicating that neither the new wind stress parameterization nor the wave-current interaction has any effect when the waves are small.
Along the other flight track sections during the 9 September flight and the entire flight on 12 and 14-15 Sep-tember, the results show that experiment A significantly overestimates H s almost everywhere (Figs. 4d, 5d, and 6d), and the error increases as H s increases (Fig. 7a). The H s prediction in experiment B is generally lower compared to experiment A. Take the 9 September flight, for example; the root-mean-square error of H s in Table 1 is reduced from 2.25 m in experiment A to 1.67 m in experiment B (the reduction is about 26%). This is because the new parameterization reduces the wind stress at higher wind speeds, and hence reduces the wind input to waves in the model. However, the H s values are still considerably larger than observations (Fig. 7b).
When the ocean current is introduced to the wave model in experiment C, the root-mean-square error is further reduced to 0.9 m for the 9 September flight (about 60% error reduction), the overall agreement with the observations is significantly improved, and the systematic overprediction for high wind speeds has been removed (Fig. 7c). This is consistent with the finding in Fan et al. (2009a). They investigated the surface wave and ocean current responses under idealized tropical cyclones, and also found that the wave-current interaction tends to significantly reduce the magnitude of simulated H s . In the next subsection, a detailed analysis is given to The H s comparison for the 9 September flight in Fig. 4d shows that the results from experiment C almost overlap with the observations everywhere except for a small section along the flight track in front of the hurricane (shown by the red dashed line in Fig. 4a). Although the H s prediction is significantly improved in this section compared to experiments A and B, the differences between the model results and observations are still significant. Furthermore, the H s comparison along the 12 September flight also showed large differences between the model results from experiment C and observations along most sections of the flight track. However, the H s comparison for the 14-15 September flight does not seem to have this problem. This discrepancy could be due to the influence of the preexisting mesoscale variability of the ocean current in the Caribbean Sea on surface waves, which is not considered in this study. Our ocean model initialization methodology provides a realistic representation of the Loop Current and mesoscale eddies in the Gulf of Mexico, but no real-time data assimilation is done in the Caribbean Sea. Instead, the GDEM monthly climatology data are used to initialize the 3D temperature and current fields. Since the climatology data smooth out most of the mesoscale features, the modeled current field also shows a smooth structure in the Caribbean area. The effect of mesoscale features, such as the Loop Current and a warm-core ring in the Gulf of Mexico, on the wave predictions will be discussed in detail in section 5d.
The model results in all three experiments indicate consistent underestimation of H s within the hurricane eye region, except on 9 September when the radius of maximum wind was small (13 km) and Hurricane Ivan was moving with a relatively fast forward speed (6 m s 21 ). On 14-15 September, when the radius of maximum wind was 3 times larger (39 km) and Ivan was moving about 3 times slower (2 m s 21 ), significant downward excursions in the model H s for each of the six times the aircraft passed through the eye are shown (Figs. 6d, 14). We speculate that a possible reason for the degraded H s predictions near the hurricane eye might be the not very accurate representation of the inflow angle in the HRD surface wind fields in these instances. From idealized experiments (not shown here), we found that variations of the inflow angle in the surface wind have a very small effect on the surface wave prediction when the hurricane has a small eye and moves relatively fast, like Hurricane Ivan on 9 September. But for slowly moving hurricanes with large eyes, the H s prediction can be significantly affected by even small changes in the wind inflow angle. HRD surface wind analysis is based on available surface wind observations from buoys, Coastal-Marine Automated Network (CMAN) platforms, ships, and other surface facilities. Because these data are often sparse near hurricanes, aircraft flight-level observations adjusted to the surface with a planetary boundary layer model (Powell 1980) are used to supplement the in situ surface measurements. Based on examination of the in-flow angle change from the launch levels to the surface, the wind directions in HWIND for surface-adjusted flight-level winds over water are given a constant angle of about 208 (Powell et al. 1996). It is possible that the real inflow angles in the Ivan surface wind field near the eyewall were quite different from the values assigned by HWIND.
The model H s from all three different experiments is also compared with NDBC buoy data from 13 to 16 September in the Gulf of Mexico (Fig. 8). Overall, the WW3 simulations show good agreement with observations. At buoys 42003 and 42039, WW3 with the original C d parameterization (experiment A) overestimates the maximum H s by about 1.5-2 m, while the simulations in experiments B and C yield much reduced errors. Buoy 42040 was adrift after 2100 UTC 15 September, which introduces some uncertainty in the accuracy of the comparison with the data. Over all, despite the buoy drift, the observations are in reasonably good agreement with the model predictions. On 15 September, Hurricane Ivan passed directly over six wave-tide gauges deployed by the Naval Research Laboratory north of buoy 42040 (;29.38N, 878W). Three buoys observed waves with H s reaching maximum values of 17.9, 16.1, and 17.1 m (Wang et al. 2005). These values are in good agreement with the model-produced H s in experiment C, as seen in the significant wave height wave swath (Fig. 10). Figure 9 shows H s measured by the radar altimeters on the Envisat-1 and ERS-2 satellites at 0338 and 0406 UTC 15 September, respectively (tandem track shown in Fig. 1), compared with WW3 results from experiment C at 2200 UTC 14 September and 0200 and 0400 UTC 15 September. The open circles show SRA H s observations between 2104 UTC 14 September and 0257 UTC 15 September, which were within 10 km of the satellite track. The simulated H s on 15 September, the closest time to the satellite measurements, compares well to the satellite data, although the altimeter shows higher maximum H s values by 1-2 m compared to both the model results and SRA data. The satellite observations and model predictions in Fig. 9 indicate spatial variations at particular times. The spatial variation of the SRA data, collected over a 6-h interval, should not be expected to match any particular model curve. For example, the four SRA data points clustered near 170 km and 6.6 m were acquired at about 2235 UTC 14 September. They are slightly above the dotted model curve for 2200 UTC. On the other hand, the three SRA data points clustered at about 150 km and 8 m were acquired at about 0220 UTC 15 September and are quite close to the dashed model curve for 0200 UTC. The comparison between the predicted H s and satellite measurements confirms that including wave-current interaction improves wave forecast skill.
The swath pictures of H s in all three experiments are shown in Fig. 10 (Figs. 4c, 5c, 6c). The same tendency was found for the H s simulations. Unlike the H s results, however, the dominant wave lengths are noticeably shorter than those in the SRA measurements. To the right of the hurricane, where the waves are more developed, the three experiments yield very different DWL values (Figs. 4c,5c,6c). In experiment A, DWL are mostly longer than those in the SRA observations (Fig.  7d), while in experiment C they are shorter (Fig. 7f).  Overall, WW3 seems to underestimate DWL when H s is correctly predicted. Fan et al. (2009b, manuscript submitted to J. Phys. Oceanogr.) have also noticed this tendency when their results are compared with DWL empirical formulas from other studies. The underestimation of DWL is most likely due to the nonlinear wave interaction term calculated within WW3. The deficiencies of the discrete interaction (DIA) are discussed in detail in Vledder et al. (2001). Tolman (2004;H. Tolman 2008, personal communication) also shows that the present WW3 nonlinear interaction calculation based on the DIA systematically overestimates the wind sea spectral peak frequency by roughly 10% (i.e., underestimates the DWL of wind seas by roughly 20%). For swell, such biases are not obvious in WW3.

b. Reduction of significant wave height by ocean currents
A significant finding in the previous section is that inclusion of the ocean current systematically reduces H s prediction. There are two ways the ocean current impacts the wave field as described in the methodology section, that is, the subtraction of the current vector from the 10-m wind vector, and the modification of the wave action equation. Both effects are included in experiment C. To determine which current effect is more important, experiment D is designed. This experiment is the same as experiment C, except the effect of current in the wave action equation is not considered. The results of the wave field simulation corresponding to the time/location of the SRA measurement at 1800 UTC 9 September are presented in Fig. 11. Even though we presented a snapshot in Fig. 11, these results are very representative throughout the whole flight period.
Starting from experiment B (without current effects), if the 10-m wind speed input is modified by the current but the wave action equation is not affected (experiment D), the resulting simulation of H s indicates small changes, as seen in Fig. 11d. Notice that H s is reduced in the area where the wind and current vectors have similar directions and increased where the wind and current vectors misaligned, as seen from Figs. 11b and 11d. When the current effect in the wave action equation is also included (experiment C), H s is significantly reduced, especially where H s reaches its maximum, as shown in Fig. 11c. These figures clearly indicate that the current effect on the wave field is mainly through the wave action equation. The relative wind speed effect is significantly smaller.
Let us next examine why including the ocean current in the wave action equation tends to reduce H s . Since the direction of the dominant wave is mostly within 308 of the direction of the current (Figs. 11, 12), we can consider for simplicity a one-dimensional approximation of the wave action equation. Furthermore, the wave action equation is expressed in the coordinate system moving with the hurricane, and the time tendency term is neglected (i.e., the wave field is assumed stationary in the moving coordinate). Then, the wave action equation used in experiment C [Eq. (1)] is simplified to where C g is group velocity and U t is the hurricane translation speed projected onto the wave propagation direction s. If we only consider the current effect on relative wind speed (experiment D), then (2) is further simplified to where the subscript 0 in N and F denotes no current.
Subtracting (2) Then, Eq. (4) shows that the reduction of the wave action spectrum (N 0 2 N) from experiment D to experiment C, shown in Fig. 11c, is caused by three factors. First, when waves are compressed or stretched by a spatially varying current, the resulting modulation of the wave action is expressed by the term 2k(›N/›k)(›U c /›s). Second, the term (›N/›s)U c is the modulation to the wave field due to horizontal current advection. This term can be interpreted as follows: if the forcing term is set such that the wave field grows with fetch (›N/›s . 0), then the spatial wave growth is reduced by a positive current simply because the wave packet propagates faster. The third effect is the modification of the forcing term (F 0 2 F), which is expected to be more important for shorter waves (spectral tail). Let us consider a wave pathway (pink arrow) in Fig. 11. Along this path the reduction of H s [i.e., (N 0 2 N) near the spectral peak] rapidly increases (Fig. 11c). Along the same path, H s (and therefore N near the spectral peak) increases (Fig. 11a) and the ocean current U c remains large (Fig. 11b). After close examination of the spectral output along this path (not shown), we have found that the compression/stretching term is relatively unimportant near the spectral peak and that the advection term (›N/›s)U c is mainly responsible for the reduction of the significant wave height (i.e., waves become lower when the wave group propagates faster because of the positive ocean current). We have also examined the magnitude of all terms in the full (2D) wave action equation, and have confirmed that the advection term along the wave propagation is dominant over a large area where the current is strong and roughly aligned with the wave propagation direction, yielding the significant reduction of H s . This analysis also highlights the significance of the hurricane translation speed U t . Equation (4) indicates that the reduction of N is enhanced as (C g 2 U t ) decreases (i.e., as the translation speed increases). In fact, when C g , which is typically ;9-10 m s 21 , is close to U t (near resonance), the reduction of N becomes the largest.

c. Wave spectrum
Next, we compare individual model spectra obtained at various positions along the 9 September flight track with those of the SRA measurements. Five spectra are selected for the comparisons (white points A-E shown in Fig. 4a). Figure 12 shows the SRA directional wave spectra and the model spectra in all three experiments at locations A-E. All three experiments show good agreement with the observations in simulating the peak wave direction. From Fig. 4a, we can see that locations A-D are in front of the hurricane, and the waves there are actually swells propagated in the tangential direction from the radius of maximum wind at an earlier position of the storm. They were generated because of the resonance, that is, they were exposed to prolonged forcing from wind because the hurricane translation speed was comparable to the group speed of the dominant waves (Moon et al. 2003;Young 2006). We can see that at all five locations, the model produces higher peak energy in experiments A and B, but similar peak energy to observations in experiment C. Also notice that the angular distribution of the wave energy in experiment C is widened. The directional spreading tends to become wider when the ocean current is included in the WW3 simulation, being consistent with the Tolman et al. (1996) study of wave interference with the Gulf Stream. This is likely caused by spatial variation of the ocean current, although it is difficult to quantitatively examine the current effect on directional spreading.
At locations A, C, and D, the model produces narrower directional spreading than in the observations in experiments A and B, but similar directional spreading to observations in experiment C. However, at locations B and E, the model produces similar directional spreading in experiment B, but larger directional spreading in experiment C. This discrepancy may be related to some overestimation of the ocean current. In this study, we used the bulk formula for calculating wind stress in POM with the drag coefficient parameterization based on the CWW model. Fan et al. (2009b, manuscript submitted to J. Phys. Oceanogr.) have pointed out that the momentum flux into the ocean can be significantly reduced because of the spatial and temporal variations of the hurricaneinduced surface waves. Fan et al. (2009a) have also shown that the coupled wind-wave-current processes can significantly reduce the momentum flux into the ocean in the right rear quadrant of the hurricane. Since these processes are not considered in our experiments, the momentum flux input to the ocean is likely overestimated. As a result, the currents and current divergence are overestimated too, especially at locations to the right of the hurricane track. Because both B and E are located close to the right of the hurricane track, the overestimation of the directional spreading in the model may be caused by the overestimation of the current. Another possibility is that, as Holthuijsen and Tolman (1991) pointed out, the existence of counter-or following current jet may modify the directional spreading of the wave spectrum. As we have discussed in section 5a, we used the GDEM monthly climatology to initialize the 3D temperature and current fields in our ocean model. Since the climatology data smooth out most of the mesoscale features, the modeled current field also shows a smooth structure in the Caribbean area and wipes out the effect of mesoscale eddies. The frequency spectra at locations A-E are shown in Fig. 13. We can see that the frequency spectrums in experiments A and B are much higher than the observations at all five locations. When the wave-current interaction is introduced in experiment C, the peak of the frequency spectrum is reduced, which greatly improves the comparison of overall (integrated) energy with observations, although it also consistently shifts the peak toward higher frequency.

d. Effect of loop current on wave prediction
To investigate the effect of preexisting currents due to mesoscale ocean features on wave prediction, we modified experiment C such that the Loop Current and its warm-core ring in the Gulf of Mexico are removed from the ocean initialization. Figure 14c shows an H s comparison between experiment C results with and without the Loop Current initialization along the 14-15 September flight. The SRA measurements are also shown for reference. The H s difference between the two experiments is clearly seen along some of the flight sections. Let us examine two such periods highlighted by the gray areas in Fig. 14c.
At 2100 UTC 14 September, H s is significantly larger with the Loop Current initialization. The spatial snapshot of the H s difference with and without the Loop Current initialization is shown at the corresponding time in Fig. 14a. Figure 15c shows the spatial distribution of the ocean temperature at 70-m depth and current field at L/(4p) depth also at the same time. At this time the aircraft is over the edge of the Loop Current (Fig. 15c), where a strong northward current is added due to the LC initialization (Fig. 16a). The wave field at the same time (Fig. 16b) indicates that the dominant waves are propagating southward at this location. If we consider the evolution history of these dominant waves (along the pink arrows in Figs. 16a,b), it is evident that a strong opposing current persisted (i.e., the packet propagation was slower) throughout the wave evolution such that the overall wave spectrum was enhanced. This explains why the predicted H s at this location is increased when the Loop Current initialization is included.
At 0240 UTC 15 September, the predicted H s is significantly smaller with the Loop Current initialization (Fig. 14c). Figure 15d shows that the flight is passing through the southern edge of the warm-core ring at this time. Because of the initialization of the warm-core ring, a strong westward current is added at that location (Fig. 16c). The wave field at the same time (Fig. 16d) shows that the dominant waves are propagating westward. Again, the evolution history of these dominant waves (along the pink arrows in Figs. 16c,d) is such that a strong positive (aligned) current accelerated the wave packet propagation and reduced the spectral level throughout the wave evolution.
These two examples clearly demonstrate that strong currents due to preexisting mesoscale ocean features may significantly modify the wave field prediction mainly because such currents accelerate or decelerate the wave propagation.

Summary and conclusions
It has been shown in previous studies that the operational wave model WW3 overestimates the significant wave height under very high wind conditions, such as under strong hurricanes. In this study we have investigated how the performance of WW3 is affected by different drag coefficient parameterizations and by including the effect of wave-current interaction. Hurricane Ivan in 2004 has been used as a test case since several observations were available for comparison, including the detailed direct observations of wave spectra from the NASA Scanning Radar Altimeter, NDBC buoy, and satellite measurements.
The drag coefficient has been parameterized by either using the original formulation in WAVEWATCH III or the coupled wave-wind model, which is based on the explicit integration of the waveform drag. The effect of wave-current interaction has been included by passing the hurricane-induced currents calculated by the Princeton Ocean Model into the coupled wave-wind model. The real-time wind analysis during Hurricane Ivan produced by the NOAA/Hurricane Research Division has been used to force both the wave model and the ocean model.
The results can be summarized as follows: 1) All experiments in this study show good prediction of wave direction, indicating that the effects of the wind stress parameterization and wave-current interaction on wave direction prediction are negligible. 2) The original WAVEWATCH III drag parameterization tends to overestimate the significant wave height, wave energy, and the dominant wavelength under very strong wind forcing, and the error seems to increase as the significant wave height increases.
3) The improved stress parameterization, together with the wave-current interaction, is shown to improve forecasts of significant wave height and wave energy. 4) The hurricane-induced ocean current tends to reduce the significant wave height mainly because it increases the advection velocity of the wave packet. Spatial variation of the current widens the directional spreading of the wave spectrum. 5) When the hurricane moves over a preexisting mesoscale ocean feature, such as the Loop Current in the Gulf of Mexico or a warm-and cold-core ring, the wave field may be significantly modified. This is mainly because strong currents associated with these features accelerate or decelerate the wave propagation and thus cause the modulation of the wave spectrum.
The results presented in this paper confirm that a fully coupled wind-wave-ocean system as suggested in Fan et al. (2009a) is necessary to accurately forecast wave fields in hurricanes.