Gulf Eddies

class: center, middle

.title[Coherent Eddy Properties from Individual Trajectories, with Application to the Gulf of Mexico] 
  
 
  
  
  
.author[Jonathan Lilly1] 
 
  
.coauthor[Paula Pérez-Brunius2, Jeffrey Early3, Julien Jouanno4, Julio Sheinbaum2]
  
.institution[1Theiss Research, 2Centro de Investigación Científica y de Educación Superior de Ensenada, 3NorthWest Research Associates, 4Laboratory for Studies in Geophysics and Spatial Oceanography] 
  
 
.date[June 18, 2019] 
  
.center[[{www.jmlilly.net}](http://www.jmlilly.net)] 
 
.footnote[Created with [{Remark.js}](http://remarkjs.com/) using [{Markdown}](https://daringfireball.net/projects/markdown/) + [{MathJax}](https://www.mathjax.org/)]

---
class: left
## Estimating Kinematic Properties from Individual Trajectories

There is a case in which this can work very well, and this is the case in which the particle is trapped inside an eddy.

The method to accomplish this draws from signal analysis.

An AM/FM signal sends information through the amplitude or frequency modulation of a carrier wave.

The trajectory of a particle trapped in a coherent eddy is a modulated oscillation in two dimensions—a stereo signal.

If the modulation is slow compared to the orbital motion, eddy properties can be recovered, as can enclosed kinematic properties.

<!--More specifically, the trajectory of a particle in the frame of an eddy is conceptualized an orbiting a time-varying ellipse.

We need to estimate the time-varying ellipse properties from a trajectory in which the eddy is superposed on other variability. -->

In the past, when datasets were small, eddy analysis involved a lot of analyst supervision and subjective.

Now we need new methods to objectively handle vast datasets, including significance tests and error quantification.

No special deployment needed, high data density not needed.

---
class: center
##Lagrangian Analysis of a Nonlinear Eddy 
 
<video preload="auto" width="95%" height="auto" data-setup="{}" autoplay loop controls><source src="../videos/nonlineareddyparticlemoviewithellipses.mp4" type="video/mp4" /></video>

An eddy at 24°N tracked by 512×256=131,072 particles, by J. Early.

The ellipses represent a time-varying best fit to local oscillatory variability in two dimensions. Color is estimated enclosed vorticity.

See the [{2018 Ocean Sciences talk}](http://www.jmlilly.net/talks/lilly18-os.html) for more details.

---
class: center
##Estimated vs. Actual Enclosed Vorticity 
 
<img style="width:70%;margin-top:-0.2em;margin-bottom:-0.2em" src="../figures/nonlineareddygoodvorticity.png">

Log10 histogram of estimated vs. true enclosed vorticity.

---
class: center
##Estimated vs. Actual Normal Strain 
 
<img style="width:70%;margin-top:-0.2em;margin-bottom:-0.2em" src="../figures/nonlineareddygoodnormalstrain.png">
 
Log10 histogram of estimated vs. true enclosed normal strain.

---
class: center
## A Drifter Dataset for the Gulf of Mexico
 
<img style="width:85%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_trajectories.jpg">

3310 drifter trajectories from 11 different sources, with 3 different drogue depths, for a total of 4.2 M hourly data points. 
---
class: center
## A Mean Flow Map
 
<img style="width:95%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_meanflowsig.jpg">

Forming a mean flow through an obvious bin average is unsatisfactory due to the temporal inhomogeneity of sampling.

---
class: center
## A Better Mean Flow Map
 
<img style="width:95%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_bettermeanflowsig.jpg">

Taking the average over each of 306 different calendar months (June 1993 through Nov 2018), then averaging again, is much better.

---
class: left
##A Model for Lagrangian Trajectories

A model for Lagrangian trajectories `$z(t)=x(t) + \mathrm{i} y(t)$` with eddies:

`\[z(t)=x(t) + \mathrm{i} y(t) =z_\epsilon(t) + \boxed{z_\star(t)}.\]`

Turbulent background `$z_\epsilon(t)$` (stochastic) plus `$z_\star(t)$` (oscillatory).         
        
The oscillatory signal `$z_\star(t)$` is modeled as a *modulated ellipse*:

`\[z_\star(t)=\mathrm{e}^{\mathrm{i}\theta(t)}\left[a(t)\cos \phi(t)+\mathrm{i}b(t)\sin \phi(t)\right].\]`

<img style="width:90%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/ellipsestability.png"> 
---
class: left
##Lagrangian Estimation of Eddy Properties

This problem has five parts. 
 
2. From `$\widehat{z}_\star(t)$`, estimate `$a(t)$`, `$b(t)$`, `$\theta(t)$`, and `$\phi(t)$` that generated it.  
2.  Address errors in `$\widehat{z}_\star(t)$` due to modulation strength, or bias error.
1.  Address “false positives” due to inertial oscillations.
3.  Address “false positives” of  `$\widehat{z}_\star(t)$` due to the background `$z_\epsilon(t)$`.
1.  Estimated the eddy signal  `$\widehat{z}_\star(t)$` given a noisy time series `$z(t)$`.

Solutions:

1. The analytic signal of .cite[Gabor (1946)], see  .cite[Lilly and Olhede (2009)].
2. The internal bias estimate of .cite[Lilly and Olhede (2012)].
1. Easy; specify $\omega/f>-1/2$, as shown shortly. 
3. A significance test based on a stochastic model, introduced here.
4. Multivariate wavelet ridge analysis of .cite[Lilly and Olhede (2012)].
    
---
class: left
## Inferring Ellipse Properties

`\[z_\star(t)=\mathrm{e}^{\mathrm{i}\theta(t)}\left[a(t)\cos \phi(t)+\mathrm{i}b(t)\sin \phi(t)\right]\]`

Given  $z_\star(t)$, we wish to infer `$a(t)$`, `$b(t)$`, `$\theta(t)$`, and `$\phi(t)$`.

This is a generalization of, given  `$x_\star(t)=a(t)\cos \phi(t)$`, infer  `$a(t)$` and `$\phi(t)$`, solved in a landmark paper by  .cite[Gabor (1946)].

The solution is to form the analytic signal `$x_\star^+(t)$`, by zeroing negative frequencies and doubling positive frequencies.

Then  `$x_\star(t)=\Re\left\{x_\star^+(t)\right\}$`, and `$\widehat{a}(t)\mathrm{e}^{\mathrm{i}\widehat{\phi}(t)}\equiv x_\star^+(t)$` gives sensible estimates.
This inverse problem has no other sensible solution .cite[(Vakman 1996)].

Generalizing to two variables, we form `$\mathbf{z}(t)=[x(t)\,\,y(t)]^T$` and take its analytic part  `$\mathbf{z}_+(t)$` .cite[(Lilly and Olhede, 2009)]. This defines `$\widehat{a}(t)$`, `$\widehat{b}(t)$`, `$\widehat{\theta}(t)$`, and `$\widehat{\phi}(t)$` which are again the only sensible estimates.

“Sensible” includes “Gives the right answer for a pure sinusoid.”
     
However, when noise is present, we can't apply this method directly. The analytic signal is implicitly employed in wavelet ridge analysis.

---
class: center
## Example of Wavelet Ridge Analysis
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_exampleeddyplot.jpg">

Modulated oscillations are extracted from a trajectory in the southwest Gulf of Mexico, the Bay of Campeche.

On the right are the original, extracted, and residual latitude signals.
 
Define `$Ro\equiv \omega/f$` as a Lagrangian frequency Rossby number.          
Define `$L_\star$` as ridge length in number of cycles completed.

Ridges are extracted in the band `$1/64<Ro <2$` and with `$L_\star>1$`.

A fourth important parameter controls the time-frequency tradeoff.  
We set the number of oscillations in the wavelet center to `$\approx 1.6$`.

---
class: left
## Example of Wavelet Ridge Analysis
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_exampletransform_fortalk.jpg">

An elliptical signal is estimated by `$\widehat{\mathbf{z}}(t)\equiv \mathbf{w}\left(t,\widehat{s}(t)\right)$` where `$\widehat{s}(t)$` are “ridges”,  maxima curves of the wavelet transform modulus `$\|\mathbf{w}(t,s)\|$`.

Constrained on physical grounds to end at `$|w_+(t,s)|=|w_-(t,s)|$`.

---
class: center
## All Ridge Properties
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_rossbyscatterplot_fortalk.jpg">

Scatter plots of ridge-averaged quantities for all ridges.    
`$R_\star$` vs. `$V_\star$` (left) and `$R_\star$` vs. `$Ro_\star$` (right).

1. Most (`$\approx 90$` %) ridges are inertial oscillations
2. Massive cyclone / anticyclone asymmetry
2. Gap in anticyclonic events around `$Ro_\star=-1/2$`
3. Unclear what is potentially due to noise

---
class: left
## Removing Inertial Oscillations

Inertial oscillations generally occur at  $\omega=-f$, but may be shifted by background vorticity to $\omega=-f-\zeta/2$ (Kunze, 1985).

The strongest possible nonlinear anticyclonic eddy is $\zeta/f=-1$  (e.g. Thomas, 2007), implying $\omega/f=-1/2$ for a circular eddy.

This most nonlinear anticyclone frequency-shifts inertial oscillations from $\omega=-f$ to $\omega=-f-\zeta/2$ or $\omega/f=-1/2$.

Thus, no eddies should occur below  $\omega/f=-1/2$, and no inertial oscillations should occur above  $\omega/f=-1/2$.

There is no comparable stability boundary on the cyclonic side.

Thus we keep only ridges with ridge-averaged Rossby numbers `\[Ro_\star\equiv\omega_\star/f_\star>-1/2.\]`
Next we turn to examining false positives emerging from the stochastic background. 
---
class: center
## A Noise Dataset & Null Hypothesis
 
<img style="width:70%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_noisespaghetti.jpg">

To assess significance, a noise dataset is created.  
Isotropic velocity spectrum `$S_{\varepsilon\varepsilon}(\omega) \equiv \min\left\{\widehat S_{++}(\omega),\widehat S_{--}(\omega) \right\}$`.   
Same variance, initial location, and mean velocity as original data.

---
class: center
## A Significance Test
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_doublesignificanceplot.jpg">

`$\zeta_\star$` is the ridge-averaged circular polarization `$\frac{2ab}{a^2+b^2}$`.   
`$L_\star$` again the ridge duration in number of cycles.

Significance is assessed using survival functions on the `$\zeta_\star/L_\star$` plane.

A 90% significance curve means more extreme (larger `$|\zeta_\star|$`or `$L_\star$`) events occur 10x more often in the data than in the noise. 
     
---
class: left
## A Bias Estimate

The leading-order term in the bias error is .cite[(Lilly and Olhede, 2012)]

`\[\frac{\|\widehat{\mathbf{z}}_+(t)-\mathbf{z}_+(t)\|}{\left\|\mathbf{z}_+(t)\right\|} 
\approx \varsigma(t)\equiv \frac{1}{2} \frac{P}{\omega^2(t)}  \frac{\left\|\mathbf{c}(t)\right\|}{\left\|\mathbf{z}_+(t)\right\|}\]`

where the parameter `$P$` controls the wavelet time/frequency tradeoff.

The second-order departure  of `$\mathbf{z}_+(t)$` from a local oscillation at the frequency $\omega(t)$ is described by
`\[
\mathbf{c}(t)  = \frac{\mathrm{d}^2}{\mathrm{d}t^2}\,\mathbf{z}_+(t)-\mathrm{i}2\omega(t) \frac{\mathrm{d}}{\mathrm{d}t}\,\mathbf{z}_+(t)-\omega^2(t)\mathbf{z}_+(t)
\]`
where the *joint instantaneous frequency* is the natural measure of time-varying frequency in the bivariate case .cite[(Lilly and Olhede, 2009)]
`\[\omega(t)\equiv\frac{\Im\left\{\mathbf{z}_+^H(t)\frac{\mathrm{d}}{\mathrm{d}t}\mathbf{z}_+(t)\right\}}{\|\mathbf{z}_+(t)\|^ 2}.\]`

This states that error is proportional to the second-order deviation of the signal from a pure sinusoid over the time scale of the wavelet.

---
class: center
## Bias Near Ridge Edges
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_edgedistributions_agezeta_fortalk.jpg">

Left is histogram of all points on edge distance `$D$` vs. circularity `$\zeta$`.
Right is mean value of estimated bias `$\varsigma$` over this plane.

Low circularities and high errors occur near edges. 
Suggests trimming `$D<1/2$` to remove wavelet spin-up effects.

---
class: center
## Eddy Ridge Properties
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_eddyscatterplot_fortalk.jpg">

Statistically significant, non-inertial, trimmed ridges.    
Scatter plots of `$R_\star$` vs. `$V_\star$` (left) and `$R_\star$` vs. `$Ro_\star$` (right).

Major cyclone / anticyclone asymmetries:
1. Highly nonlinear (submesoscale?) cyclones for `$R_\star<10$` km
2. Nonlinear cyclones completely dominate for `$10<R_\star<50$` km 
3. Anticyclones (Loop Current Eddies) dominate for `$R_\star > 50$`

---
class: center
## All Anticyclones Colored by `$Ro_\star$`
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_allanticycloneswithrossbycolor.jpg">

Mostly we are seeing the well-known Loop Current Eddies.

---
class: center
## All Cyclones Colored by `$Ro_\star$`
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_allcycloneswithrossbycolor.jpg">

Many medium and small high `$Ro_\star$` events, but few large events.

---
class: center
## All Anticyclones Colored by `$L_\star$`
 
<img style="width:99%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_allanticycloneswithlengthcolor.jpg">

Anticyclones are typically short to medium duration.

---
class: center
## All Cyclones Colored by `$L_\star$`
 
<img style="width:99%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_allcycloneswithlengthcolor.jpg">

High `$Ro_\star$`, medium-sized cyclones are typically long duration.

---
class: center
## Case Study of a Loop Current Eddy
 
<img style="width:80%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_anticyclones2013.jpg">

The origin and evolution of the Loop Current Eddies is well known.   
Here we see evidence of three different propagation velocities.

---
class: center
## Strong Cyclones `$|Ro_\star|>1/6$` Colored by `$Ro_\star$`
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_strongcycloneswithrossbycolor.jpg">

The most nonlinear cyclones are in the east, west, and south.

---
class: center
## Strong Cyclones `$|Ro_\star|>1/6$` Colored by `$L_\star$`
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_strongcycloneswithlengthcolor.jpg"> 
 
Prominent medium-sized, long-lived cyclones. Are these related?

---
class: center
## Strong Anticyclones `$|Ro_\star|>1/6$` 
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_stronganticycloneswithrossbycolor.jpg">

Anticyclonic events are almost entirely absent in this band.

---
class: center
## Origin of Strong Mesoscale Cyclones
 
<video preload="auto" width="80%" height="auto" data-setup="{}" autoplay loop controls><source src="../videos/origincyclonemovie.mp4" type="video/mp4" /></video>

Lagrangian ellipses combined with sea surface height from AVISO.

---
class: center
## Strong Cyclones Feed the Campeche Gyre
 
<video preload="auto" width="80%" height="auto" data-setup="{}" autoplay loop controls><source src="../videos/mergercyclonemovie.mp4" type="video/mp4" /></video>

Lagrangian ellipses combined with sea surface height from AVISO.

---
class: center
## Comparison with a Numerical Model
 
<video preload="auto" width="72%" height="auto" data-setup="{}" autoplay loop controls><source src="../videos/roms_vorticitymovie.mp4" type="video/mp4" /></video>

Source of positive vorticity is inflow through the Yucatán Channel.

---
class: center
## Small Cyclone Ellipses, `$R_\star<10$` km
 
<img style="width:100%;margin-top:-0.7em;margin-bottom:-0.8em" src="../figures/sener_smallcycloneswithlengthcolor.jpg">

Focus on small cyclones, with 4x as many ellipses as before.

---
class: center
## Small Cyclone Trajectories, `$R_\star<10$` km
 
<img style="width:90%;margin-top:-0.7em;margin-bottom:-0.8em;margin-left:-3em" src="../figures/sener_smallcyclonesoriginal.jpg">

Trajectories often cusp rather than loop.  Large-scale advection.

---
class: center
## Small Cyclone Residuals, `$R_\star<10$` km
 
<img style="width:90%;margin-top:-0.7em;margin-bottom:-0.8em;margin-left:-3em" src="../figures/sener_smallcyclonesresidual.jpg">

Smooth residuals indicate signals, while small, are well captured. 
  
---
class: left        
###Results

An objective analysis of modulated oscillations in drifter trajectories in the Gulf of Mexico reveals new features of the eddy field.

1. Intense cyclones with `$Ro_\star>1/6$` throughout the Gulf, with no analogue on the anticyclonic side
2. Strong `$10<R_\star<50$` km cyclones in the east are Loop Current Frontal Eddies (LCFE's), e.g. Le Hénaff et al. (2014)
3. Strong `$10<R_\star<50$` km cyclones in the west appear to be LCFE's that have propagated across the basin
4. Western cyclones feed the Campeche Gyre, suggesting a vorticity origin in the cyclonic shear zone of Loop Current
5. Highly nonlinear, very small `$R_\star<10$` cyclones appear to indicate an active, rapidly propagating submesoscale 
6. Substantial anticyclone activity is not seen until the Loop Current Eddies at `$R_\star>50$`

To do: compare with eddies in numerical models of the Gulf, detected using this method and also standard Eulerian methods

.center[Visit [{www.jmlilly.net}](http://www.jmlilly.net)  for papers, this talk, and  Matlab code.]

.center[.footnote[P.S. Like the way this presentation looks? Check out [{Liminal}](https://github.com/jonathanlilly/liminal).]]