\section{Introduction} \label{sec:introduction}

Since its discovery more than a century ago~\cite{Hertz1887}, the photoelectric effect, i.e., the emission of an electron from an atom, molecule or extended target, due to the absorption of ionizing radiation, has played a fundamental role in our understanding of charge-transfer processes in matter. Photoionization often involves the excitation of transiently bound electronic states, which decay by emitting an electron on a time scale as small as a few femtoseconds~\cite{Wickenhauser2005,Gruson2016}. Until the end of the XX century, photoelectron spectroscopies were mostly limited to study photoemission processes in the stationary regime~\cite{Schmidt1992}. During the last two decades, however, the development of novel sources of sub-femtosecond extreme ultraviolet (XUV) light pulses has opened the way to the time-resolved study of electronic dynamics in atoms and molecules at its natural timescale~\cite{Hentschel2001,Paul2001,Sansone2006,Krausz2009}. Attosecond spectroscopy has allowed the measurements and theoretical computations of the photoionization time delay from different electronic valence shells~\cite{Schultze2010a,Kheifets10,Isinger17},
across autoionizing resonances~\cite{Argenti2010,Argenti2017,Galan16,Gruson2016,Cirelli18,douguet18}, or as a function of the emission angle~\cite{Comby16,Beaulieu16,Beaulieu17,Serov16,Cirelli18}.

Reconstruction of Attosecond Beating by Interference of Two-Photon Transitions (RABBITT) \cite{Paul2001, Muller2002, Agostini2004} is a popular technique in which an XUV attosecond-pulse train (APT) is used as a pump, in conjunction with a weak infrared (IR) probe pulse, to ionize a target, as a function of the pump-probe delay $\tau$. The APT is obtained through the process of High-Harmonic Generation (HHG) \cite{LHuillier1993,Corkum1993}, which takes place when a strong replica of the IR probe pulse is focused on an active medium.

The XUV pulse trains obtained with the HHG processes described above are generally weak (less than $10^{10}$~W/cm$^{2}$ on focus), and hence they can only be used to promote the absorption of a single harmonic XUV-photon by a target. When used in isolation, a pump APT ionizes the target giving rise to a photoelectron spectrum featuring sharp peaks in correspondence with the harmonics of the fundamental IR frequency. If the ionization takes place in the presence of a weak probe replica of the initial IR pulse, with a controllable delay with respect to the XUV train, however, two-photon ionization paths become possible \cite{Paul2001}. Besides simple one-photon absorption, the other most relevant ionization paths entail  also the exchange (emission or absorption) of one IR-photon.

In traditional RABBITT schemes, the APT is linearly polarized and its spectrum comprises only odd multiples of the fundamental IR frequency, $\omega_{\scriptscriptstyle\textsc{IR}}$. The absorption of an XUV photon from the $2n-1$ harmonics, followed by the absorption of one IR photon, and the absorption of an XUV photon from the $2n+1$ harmonics, followed by the stimulated emission of one IR photon, result in transition amplitudes to the same final energy in the continuum, $\mathcal{A}^{(+)}_{2n-1}$ and $\mathcal{A}^{(-)}_{2n+1}$, respectively, which give rise to the $2n$-sideband signal in the photoelectron spectrum, as  schematically illustrated in Fig.~\ref{fig:1}a. Thanks to the inherent coherence between the XUV and IR pulses, these two amplitudes interfere. Each exchange of an IR photon imparts to the ionization amplitude a phase factor $\exp(\mp i\omega_{\scriptscriptstyle\textsc{IR}} \tau)$, depending on whether the photon is absorbed or emitted. As a consequence, the photoionization probability to the sideband $I_{2n}$ oscillates as a function of the pump-probe time delay $\tau$ at twice the frequency of the IR,  $I_{2n}\propto|\mathcal{A}^{(+)}_{2n-1}\,\mathcal{A}^{(-)}_{2n+1}|\cos(2\omega_{\scriptscriptstyle\textsc{IR}}\tau+\phi_{2n})$~\cite{Veniard1996,Jimenez2014}, where $\phi_{2n}$ is a characteristic phase shift.
The phase shift $\phi_{2n}$ incorporates both the relative phase of two consecutive harmonics  and the additional phase imparted to the photoelectron by the two-photon transition itself. In the special case in which one of the two harmonics excites a metastable state $|a\rangle$, therefore, the other harmonic can serve as a holographic reference to measure the phase of the resonant two-photon transition~\cite{Argenti2010,Argenti2017,Galan16,Gruson2016,douguet18}. Indeed, the phase shift in the sideband beating exhibits the rapid excursion characteristic of resonant amplitudes as a function of the detuning of the first harmonic from the intermediate resonance. From the phase profile it is possible to reconstruct the fast decay of the metastable state to the continuum, resolved in time~\cite{Gruson2016}.

The ability of this pioneering approach to study the time evolution of autoionizing states, however, is limited. In particular, the reconstruction of one-photon resonance decay is highly distorted because the resonance excitation amplitude is buried in a two-photon transition where it is inextricably intertwined with the contribution of virtual excitations to states across a wide range of energies. Furthermore, one-photon amplitudes are imaged to sidebands through structured continuum-continuum couplings that distort the amplitude phase as a curved mirror distorts the image of a visitor of an amusement park~\cite{Argenti2017}. Finally, the duration of IR probe pulses in attosecond experiments is typically of the order of tens of femtoseconds, i.e., comparable to or even shorter than the lifetime of most autoionizing states. As a result, the energy resolution of the already distorted resonant profile is limited by the energy resolution of the probe photon.

In the RABBITT setup illustrated above, when applied to spherically symmetric targets such as closed-shell atoms, the portions of the state function responsible for the harmonic and the sideband photoelectron amplitudes have well defined parity, i.e., odd and even, respectively, corresponding to the parity of the number of exchanged photons. As a consequence, the photoelectron distributions are symmetric upon reflection on the plane perpendicular to the polarization axis of the external pulses.
\begin{figure}[t]
\includegraphics[width=\columnwidth]{fig-1-scheme.eps}
\caption{\label{fig:1} RABBITT schemes used to study the dynamics of an autoionizing resonant state $|a\rangle$. (a) traditional RABBITT. (b) One-Two RABBITT for ionization paths interfering at the resonance, or (c) at the lower harmonic with the resonance used as an intermediate stepping stone (see text for details).}
\end{figure}
In the present work we propose an alternative interferometric scheme, which we call \emph{1-2 RABBITT},  that overcomes and quantifies the distortions inherent to the reconstruction of resonance decay with traditional RABBITT. In 1-2 RABBITT, the XUV spectrum comprises isolated pairs of linearly polarized consecutive harmonics of the probe frequency, as illustrated in Fig.~\ref{fig:1}.b. This scheme can be realized, for example, by employing an attosecond pulse train with two largely dominant consecutive odd harmonics, together with the second harmonic of the fundamental IR as a probe. Conversely, it is well known that if a small second-harmonic component is added to the strong IR pulse used in the HHG~\cite{Laurent2012, Mauritsson2008}, the resulting pulse train can comprise both even and odd multiples of the fundamental IR frequency. In the latter case, two consecutive harmonics could conceivably be isolated with the help of a suitable metallic filter. Regardless of the specific generation mechanism, for the sake of being definite, as well as to maintain a close parallel with traditional RABBITT schemes, in the following we will label the two consecutive harmonics as $2n$ and $2n+1$.
The one-photon ionization amplitude $\mathcal{A}_{2n+1}$, tuned to an isolated resonance $|a\rangle$, interferes with the two-photon amplitude $\mathcal{A}^{(+)}_{2n}$, giving rise to an asymmetric photoelectron emission along the laser polarization $\hat{\epsilon}$, $\Sigma(E)=\int d\Omega\, I_{2n}(E,\hat{\Omega})\sgn(\hat{\Omega}\cdot\hat{\epsilon})$, that beats at the IR frequency as $\Sigma(E)\propto|\mathcal{A}_{2n+1}\,\mathcal{A}^{(+)}_{2n}|\cos(\omega_{\scriptscriptstyle\textsc{IR}}\tau+\phi_{2n+1})$, where $\phi_{2n+1}=\arg\mathcal{A}_{2n+1} - \arg\mathcal{A}^{(+)}_{2n}$. In contrast to the traditional RABBITT scheme, therefore, in this case the phase of the one-photon transition to the resonance $|a\rangle$ is directly encoded in $\phi_{2n+1}$, and can be extracted from it, provided that the two-photon amplitude is not structured. Furthermore, the same experiment quantifies the distortion of the resonant profile introduced by the probe photon. Indeed, at the harmonics $2n$, the amplitudes $\mathcal{A}_{2n}$ and $\mathcal{A}^{(-)}_{2n+1}$ interfere (see panel (c) of Fig.~\ref{fig:1}),
leading to a $\omega_{\scriptscriptstyle\textsc{IR}}$-beating of $\Sigma(E)$ with terms of the form $\propto|\mathcal{A}_{2n}\,\mathcal{A}^{(-)}_{2n+1}|\cos(\omega_{\scriptscriptstyle\textsc{IR}}\tau+\phi_{2n})$, where $\phi_{2n}$ now exhibits the excursion of the intermediate resonant phase modulated by a free-free transition induced by the IR field.

In either case, the beating distribution has odd parity, since it originates from the interference between odd and even partial waves, and hence it can only be observed with a directional detector. On the other hand, the beating amplitude is proportional to the XUV intensity and to the peak IR field, which makes this scheme more sensitive than traditional RABBITT, where the signal is proportional to the IR intensity instead. Notice that the relative phase of the two consecutive harmonics is not relevant here, since the observable of interest is the rapid variation of the photoionization phase across the energy of a single harmonic, rather than its absolute value.

To illustrate the 1-2 RABBITT scheme we will employ a one-dimensional model that exhibits the key features of asymmetric resonant photoemission. The choice of a model whose attributes are analytically known allows us to better highlight the correspondence between scattering phase, photoionization phase, and the phase reconstructed from one-photon and two-photon resonant transitions, than with a many-body systems, in which the decay is driven by correlation. Furthermore, in the context of perturbative single-ionization processes, a 1D model is justified as it already reproduces most of the features of atomic ionization. This circumstance is to be contrasted with the case of strong-field ionization, where 1D models have a more limited validity since they cannot reproduce the transverse spreading of the electronic wavepacket, which is an essential aspect of all processes that depend on the recollision of the photoelectron with the parent ion.

The paper is organized as follows. In \hbox{Sec. \ref{sec:time-evolution}} we summarize the time-dependent formalism appropriate for 1-2 RABBITT, illustrate the one-dimensional model, and describe its relevant observables.
 In \hbox{Sec. \ref{sec:Results}}, we present the results for selected one-photon and two-photon transition amplitudes obtained, for several pump-probe time delays, using time-dependent perturbation theory (PT) as well as by integrating the time-dependent Schr\"{o}dinger equation (TDSE) numerically. The excellent agreement between these two approaches confirms that the analytical perturbative formulas capture all the aspects of the 1-2 RABBITT interferometric scheme. From the photoionization asymmetry computed with either methods, we retrieve the one-photon and the two-photon complex resonant ionization amplitude. In section {\ref{sec:CONCLUSION} we present our conclusions and perspectives.

\section{Theoretical approach} \label{sec:time-evolution}
\subsection{Time evolution}
Atomic units ($\hbar=1$, $e=1$, $m_e=1$) are used throughout, unless stated otherwise. The time evolution of the wavefunction $\Psi(t)$ of a quantum system driven by external time-dependent fields is governed by the TDSE,
\begin{eqnarray}
	i\frac{\partial \Psi(t)}{\partial t}=H(t)\Psi(t) \label{eq:TDSE},\\
	H(t)=H_0+H'(t) \label{eq:ham}
\end{eqnarray}
where the total Hamiltonian $H(t)$ has a field-free component, $H_0$, and a radiation-matter interaction component, $H'(t)$. In the 1D model examined in this work, $H_0=p^2/2+V(x)$, where $p=-i d/dx$ is the electron momentum and $V(x)$ is a local potential. In dipole approximation, the minimal-coupling interaction term is $H'(t)=F(t)o$, where $F(t)$ is the time-dependent field, and $o$ a suitable dipole operator. In velocity gauge $H'(t)=A(t)p$, where $A(t)$ is the vector potential, whereas in length gauge $H'(t)=-\mathcal{E}(t)x$, where $\mathcal{E}=-\partial_t A$ is the external electric field~\cite{Joachain}.
In the case of weak fields with finite duration, the solution of Eq. \ref{eq:TDSE} can be found as a truncated perturbative series. To express the perturbative solution, it is convenient to reformulate the TDSE in the interaction representation and in integral form,
\begin{equation}
	\Psi_I(t)=e^{iH_0t}\Psi(t),\qquad H^\prime_I(t) =e^{iH_0t} H^\prime(t) e^{-iH_0t},
\end{equation}\vspace{-15pt}
\begin{equation}
	\Psi_I(t)=\Psi_I(t_0)-i\int_{t_0}^t dt'H_I'(t')\Psi_I(t').
\end{equation}
The wave function is then the sum of finite-order terms,
\begin{eqnarray}
	\Psi_I(t)&=& \sum_m\Psi^{(m)}_I(t),\\
	\Psi_I^{(m)}(t)&=& -i\int_{t_0}^t dt'H_I'(t')\Psi^{(m-1)}_I(t').
\end{eqnarray}
Let us assume that the field-free Hamiltonian has a non-degenerate ground state, $H_0|g\rangle=E_g|g\rangle$, normalized to unity, $\langle g | g \rangle = 1$, and a set of ionization channels $H_0|\alpha E\rangle=E|\alpha E\rangle$, normalized as $\langle \alpha E | \beta E'\rangle=\delta_{\alpha\beta}\delta(E-E')$, where $E$ is the total energy and $\alpha$ is a set of additional quantum numbers needed to resolve possible degeneracies. If the system is initially in the ground state,  the ionization probability amplitude $\mathcal{A}_{\alpha E\leftarrow g}$ to the final state $\ket{\psi_{\alpha E}}$ at any time $t$ after the end of the pulse is
\begin{equation}
\mathcal{A}_{\alpha E\leftarrow g }=\sum_{m}\mathcal{A}^{(m)}_{\alpha E\leftarrow g}, \quad \mathcal{A}^{(m)}_{\alpha E\leftarrow g}=\braket{\psi_{\alpha E}|\Psi_I^{(m)}(t)}.
\end{equation}
The first- and second-order ionization amplitudes can be expressed in the frequency domain~\cite{Jimenez2014} as
\begin{eqnarray}
\mathcal{A}_{\alpha E\leftarrow g}^{(1)}&=&-i\braket{\psi_{\alpha E}|o|\psi_g}\tilde{F}(\omega_{Eg}), \label{eq:one-photon}\\
\mathcal{A}_{\alpha E\leftarrow g}^{(2)}&=&-i\int_{-\infty}^{\infty}d\omega\tilde{F}(\omega_{Eg}-\omega)\tilde{F}(\omega)\mathcal{M}_{\alpha E\leftarrow g}^{(2)}(\omega)\label{eq:two-photon},
\end{eqnarray}
where the two-photon matrix element has the form
\begin{eqnarray}
\mathcal{M}_{\alpha E\leftarrow g}^{(2)}(\omega)=\braket{\psi_{\alpha E}|\,o\, G_0^{+} (E_g+\omega) \,o\,|\psi_g}\label{eq:two-photon-matrix},
\end{eqnarray}
$\tilde{F}(\omega)=1/\sqrt{2\pi}\int dt F(t)\exp(-i\omega t)$ is the Fourier transform of the field, $G_0^+(E)=(E-H_0+i0^+)^{-1}$ is the retarded resolvent of the field-free Hamiltonian, and $\omega_{Eg}=E-E_g$ is the overall excitation energy~\cite{Joachain}. If the contribution of intermediate bound states is negligible (as it is certainly the case if the potential supports only one bound state), to the two-photon transition is
\begin{eqnarray}
\mathcal{M}_{\alpha E\leftarrow g}^{(2)}(\omega)=\sum_\beta\hspace{-2pt}\int\hspace{-2pt}\frac{\langle \psi_{\alpha E}|o|\psi_{\beta E'}\rangle\langle \psi_{\beta E'}|o|\psi_{g}\rangle}{E_g+\omega-E'+i0^+}dE'.\quad\label{eq:two-photon-matrix}
\end{eqnarray}
From the bound-free $\langle \psi_{\alpha E}|o|\psi_g\rangle$ and free-free $\langle \psi_{\alpha E}|o|\psi_{\beta E'}\rangle$ dipole matrix elements, therefore, it is possible to determine the photoionization probability and asymmetry resulting from the interaction between the system and an arbitrary sequence of weak pulses, under the assumption that all terms beyond second order are negligible. An alternative way to compute the ionization amplitude, which is valuable as a proof of principle to validate the perturbative results, is to solve the TDSE numerically, as discussed in the next sections.

\subsection{Analytical 1D model}
We consider a symmetric analytical model (AM) in which the potential $V(x)$ is given by two repulsive delta functions, located at $x=\pm a$, and an attractive delta potential located at the origin.
\begin{equation}
V(x)=V_+\delta(x+a)+V_+\delta(x-a)-V_-\delta(x), \label{eq:pot_analy}
\end{equation}
Where $V_+$ and $V_-$ are positive parameters that express the strength of the potential.
\begin{figure}[t]
	\includegraphics[width=0.7\columnwidth]{fig-2-states.eps}
	\caption{Analytical potential used in this study. Also shown the absolute value squared (a) of the ground state and (b) of a resonant state with odd parity.}
	\label{fig:2}
\end{figure}
\paragraph{Eigenstates}
The stationary states of the field-free Hamiltonian can be found analytically. Where the potential is zero the Hamiltonian eigenfunctions are
\begin{eqnarray}
	\psi(x)&=&Ae^{\kappa x} + Be^{-\kappa x},\quad E<0,\quad \kappa =\sqrt{-2E},\\
	\psi(x)&=&Ae^{ikx} + Be^{-ikx},\quad E>0,\quad k =\sqrt{2E}.
\end{eqnarray}
The solution must be continuous, and its derivative satisfy the well known conditions
\begin{eqnarray}\label{eq:boundary-conditions}
\psi'(0^+)&=&\psi'(0^-)-2V_-\psi(0),\\
\psi'(\pm a^+)&=&\psi'(\pm a^-)\,+\, 2\,V_+\psi(\pm a),
\end{eqnarray}
where $\pm a^{+}$ and $\pm a^{-}$ denote the right and left limit towards $x=\pm a$, respectively. Since the potential is symmetric about the origin, the parity operator $\Pi$, defined by $\Pi\ket{x}=\ket{-x}$, commutes with the field-free Hamiltonian $\left[H_0,\Pi\right]=0$. As a consequence, the eigenfunctions of $H_0$ can be assumed to be either even, $\psi_{eE}$, or odd, $\psi_{oE}$, under inversion of the spatial coordinate.

It is easy to ascertain that the potential~\eqref{eq:pot_analy} supports only one bound state.
If the repulsive potentials are far enough from the origin or, conversely, the attractive potential is sufficiently strong,
the wavefunction of the bound state is negligible at $|x|\ge a$ and so is the discontinuity of its derivative.
For a whole range of parameters, therefore, the ground state is well approximated by that of a single attractive delta potential,
\begin{equation}
\psi_g(x)=A_ge^{{-\kappa |x|}},\quad E_g=-V_-^2/2,
\end{equation}
where $A_g=\sqrt{V_-}$ and $\kappa=V_-$. In particular, the ground state is an even function, as shown in Fig.~\ref{fig:2}a.
The even scattering functions are
\begin{equation}\label{eq:phi_even}
\hspace{-3pt}\psi_{eE}(x)=\frac{1}{\sqrt{\pi k}}\left\{
\begin{array}{ll}
A_e(k)\sin[k|x|+\delta_e(k)]&|x| < a\\\vspace{-8pt}\\
\sin[k|x|+\eta_e(k)]&|x|>a
\end{array}
\right.
\end{equation}
where $\delta_e(k)=-\atan{\left(k/V_-\right)}$ is the wave-function phase inside the barrier ($|x|<a$),
$\eta_{e}(k)$ is the scattering phaseshift, and $A_e(k)$ the asymptotic amplitude,
\begin{eqnarray}
\eta_e(k)&=&\acot\left[\cot(ka+\delta_e)+2V_+/k\right]-ka\nonumber\\
        A_e(k)&=&\big[1+4V_+^2\,\sin^2[ka+\delta_e(k)]/k^2+\nonumber\\
        &&+2V_+\,\sin[2ka+2\delta_e(k)]/k\big]^{-1/2}.\nonumber
\end{eqnarray}
The odd scattering functions are
\begin{equation}\label{eq:phi_odd}
\hspace{-3pt}\psi_{o E}(x)=\frac{1}{\sqrt{\pi k}}\left\{
\begin{array}{ll}
A_o(k)\sin(kx)&|x| < a\\\vspace{-8pt}\\
\sin[kx+\sgn(x)\,\eta_o(k)]&|x|>a
\end{array}
\right.
\end{equation}
whose scattering phase and asymptotic amplitudes are
\begin{eqnarray}
\eta_o(k)&=&\acot\left[\cot(ka)+2V_+/k\right]-ka,\nonumber\\
        A_o(k)&=&\big[1+4V_+^2\sin^2(ka)/k^2+2V_+\sin(2ka)/k\big]^{-1/2}.\nonumber
\end{eqnarray}
Notice that the scattering states in~(\ref{eq:phi_even},\ref{eq:phi_odd}) are purely real, $\psi_{e/o\,E}^*(x)=\psi_{e/o\,E}(x)$, and normalized as
\begin{equation}\label{eq:normalization}
\braket{\psi_{\alpha E}|\psi_{\beta E'}}=\delta\left(E-E'\right)\delta_{\alpha\beta},\quad \alpha,\beta\in\{e,o\}.
\end{equation}
The repulsive potentials have a dramatic effect on the continuum as they can keep the electron outside, or confine it within the barrier. As an example, Fig.~\ref{fig:2}.b shows the probability density of an odd scattering state at an energy close to that of the second odd bound state of a $2a$-wide box. Since in this case the potential barriers are only partially reflective, the bound state manifests itself as a resonance embedded in the continuum. The peak density within the region $[-a,a]$ is much larger than for $|x|>a$, which shows that the potential can transiently bind an electron. The interference mechanism that gives rise to this resonance confinement is the same as in the Fabry-P\'erot etalon of classical optics~\cite{Perot1899}.
It is useful to introduce also a second set of scattering states, $|\psi_{LE}^-\rangle$ and $|\psi_{RE}^-\rangle$, which satisfy incoming boundary conditions and represent a mono-energetic electron that, \emph{after} colliding with the potential (hence the minus), emerges on the left and on the right of the potential, respectively,
\begin{eqnarray}
	\psi_{LE}^-(x)&=&\frac{1}{\sqrt{\pi k}}\left\{
	\begin{array}{ll}
		e^{-ikx}+\beta e^{+ikx} \quad& x \le -a\\
		\gamma e^{-ikx} &  x \ge a
	\end{array}\right.\\
	\psi_{RE}^-(x)&=&\frac{1}{\sqrt{\pi k}}\left\{
	\begin{array}{ll}
		\gamma e^{ikx}& x \le -a\\
		\beta e^{-ikx}+e^{ikx}  \quad& x \ge a
	\end{array}\right.
\end{eqnarray}
where the complex numbers $\beta\equiv-\left(e^{-2i\eta_e}+e^{-2i\eta_o}\right)/2$ and $\gamma \equiv -\left(e^{-2i\eta_e}-e^{-2i\eta_o}\right)/2$ satisfy $|\beta|^2+|\gamma|^2=1$.
The left and right states, which turn into one another under the action of parity, $\Pi|\psi_{LE}^-\rangle=|\psi_{RE}^-\rangle$, can be expressed in terms of the even and odd scattering states,
 \begin{equation}\label{eq:LRvseo}
       \psi_{L/R\,E}^-(x)=ie^{-i\eta_e}\psi_{eE}(x)\mp ie^{-i\eta_o}\psi_{oE}(x).
\end{equation}

\paragraph{Dipole transition matrix elements}
To evaluate the interaction of a charged particle in the model potential with an external time-dependent electric field, we must compute the dipole matrix for both bound-free and free-free transitions.
The only bound-free transition in the present 1D model is the one between the even ground state and the odd scattering states. In velocity gauge,
\begin{equation}
\mu_{E,g}=\braket{\psi_{oE}|p|\psi_g}=-i\int_{-\infty}^{\infty}\hspace{-5pt}dx \,\,\psi_{oE}(x)\frac{d}{dx}\psi_g(x).
\end{equation}
Under the assumption that the ground state wavefunction is negligible beyond $|x|>a$, the dipole integral can be evaluated using the analytical expression for both the initial and final wavefunctions in the inner region only,
\begin{eqnarray}
	\mu_{E,g}&=& -\frac{2iA_o(k)A_g\kappa k}{\sqrt{\pi k}(\kappa^2+k^2)}.
	 \label{eq:one-photon-mx}
\end{eqnarray}
Since neither the bra nor the ket are normalizable, the free-free dipole matrix element comprises both a delta and an \emph{off-shell} component,
\begin{eqnarray}
	&&\braket{\psi_{eE} |\frac{d}{dx}|\psi_{oE'}}=
	\frac{\mathcal{P}}{E'-E}\,\braket{\psi_{eE}|\left[d/dx,H_0\right]|\psi_{oE'}}\nonumber\\
	&&+\delta(k-k')\,\pi\,\sin\left\lbrack\eta_e(k)-\eta_o(k')\right\rbrack A_e(k)A_o(k'),\label{eq:ccdme}
\end{eqnarray}
where $\mathcal{P}$ indicates the principal part, and we have used the fact that the bra and ket are eigenfunctions of the field-free Hamiltonian to factor out the matrix element of $\left[d/dx,H_0\right]$, which is regular. The delta component is computed from the asymptotic expression of the scattering states. The regular factor in the \emph{off-shell} component can be readily shown to have the following expression,
\begin{equation}
\label{eq:td-cc}
\begin{split}
&\braket{\psi_{eE}|\left[d/dx,H_0\right]|\psi_{oE'}}=A_e(k)A_o(k') k k^\prime\times\\
    &\times \left\lbrace
    \cos\left[ka+\delta_e(k)\right]\cos(k' a)
    -\cos\left[\delta_e(k)\right]\right\rbrace-\\
    &-k k^\prime\cos\left[ka+\eta_e(k)\right]\cos\left[k^\prime a+\eta_o(k')\right].
\end{split}
\end{equation}
The bound-free and free-free dipole matrices can be used to compute the one- and two-photon transition matrix elements. These transition matrices are used in turn to compute the finite-pulse multi-photon integrals to determine the one- and two-photon perturbative ionization amplitudes~\eqref{eq:one-photon} and~\eqref{eq:two-photon}.

\subsection{Numerical 1D model}
In the numerical version of the model, both the potential and the wavefunction are tabulated on a grid. Whereas the singular potential $V(x)$ of Eq.~\eqref{eq:pot_analy} cannot be represented in this context, it can be approximated by a potential $V_N(x)$ in which the delta distributions are replaced with narrow Lorentzian functions $L_\Delta(x)$,
\begin{equation}
	V_N(x)=V^N_+L_\Delta(x+a)+V^N_+L_\Delta(x-a)-V^N_-L_\Delta(x)\label{eq:pot_num},
\end{equation}
where $L_\Delta(x)=\Delta/[2\pi(x^2+\Delta^2/4)]$, with $\Delta\ll a$.
The bound $\psi_g$ and scattering states $\psi_{e/o E}$ of the eigenvalue problem $H_0\psi = E\psi$ with the potential~\eqref{eq:pot_num} are computed  numerically using the renormalized Numerov method \cite{Johnson1977,Johnson1978}.
The scattering functions outside the potential barrier have the same form as the one given in~\eqref{eq:phi_even} and \eqref{eq:phi_odd} for $|x|>a$.
The constant parameters $V^N_+$ and $V^N_-$, however, may not coincide with $V_+$ and $V_-$ if the energy difference between the bound and resonant states in the NM is to match that in the AM. The values of $\eta_{e/o}(k)$ are computed numerically by matching the eigenfunction of $H_0$ to the form \eqref{eq:phi_even} and \eqref{eq:phi_odd} at a point $x_M$ outside the potential barriers.

\paragraph{Dipole transition matrix elements}
In the NM, the bound-free transition matrix elements are computed by numerical integration. For the free-free dipole matrix elements, the delta components \emph{on shell} is dictated exclusively by the asymptotic behavior of the continuum wave functions and hence it has the same form as the \emph{on shell} component of the free-free matrix element in the AM in~\eqref{eq:ccdme}.
The \emph{off-shell} component of the free-free dipole matrix element in the NM is computed by partitioning the integral $\int^{\infty}_{-\infty}\psi_{eE}(x)\psi'_{oE'}(x)dx$ into a short-range ($|x|<x_M$) and a long-range part, where $x_M$ is chosen so that the potential barrier is negligible for $|x|>x_M$. The short-range part, $2\int^{x_M}_{0}\psi_{eE}(x)\psi'_{oE'}(x)dx$, is evaluated numerically,
whereas the outer component can be expressed as the following boundary term
\begin{eqnarray}
&&\int^{\infty}_{x_M}\hspace{-3pt}\psi_{e E}\frac{d\psi_{oE'}}{dx}\,dx=\frac{k'}{2}\left(E-E^\prime\right)\times\nonumber\\
&\times\big\lbrace& k\cos\left[kx_M+\eta_e(k)\right]\cos\left[k^\prime x_M+\eta_o(k^\prime)\right]\nonumber\\
&+&k^\prime\sin\left[kx_M+\eta_e(k)\right]\sin\left[k^\prime x_M+\eta_o(k^\prime)\right]\big\rbrace.\label{eq:boundaryCorrection}
\end{eqnarray}

In the NM, the TDSE is solved in the length gauge and using finite differences with a three-point formula to evaluate the kinetic energy term. The TDSE solution is propagated forward in time using the Crank-Nicolson method~\cite{Ilchen2017}. The ionization amplitude is eventually extracted from the wavefunction by projecting it on the scattering states of the system, as $\mathcal{A}_{\alpha E\leftarrow g}=\langle \psi_{\alpha E}|\Psi(t)\rangle$.

\subsection{Physical observables}
The key observable in the 1-2 RABBITT scheme is the photoemission asymmetry, $\Sigma(E)$, defined here as the difference between the probability of photoemission to the left and to the right, per unit of energy,
\begin{equation}
\Sigma(E) = \mathcal{P}_L(E) - \mathcal{P}_R(E).
\label{eq:asymmetry}
\end{equation}
The directional photoelectron distributions are the square of the corresponding photoemission amplitudes, $\mathcal{P}_{L/R}(E)=|\mathcal{A}_{L/R\,E\gets g}|^2$, where $\mathcal{A}_{L/R\,E\gets g}=\braket{\psi^-_{L/R\,E}|\Psi_I(t)}$ are obtained by projecting $\Psi(t)$, after the end of the pulse. In turn, using~\eqref{eq:LRvseo}, the photoemission amplitudes to the left and to the right can be expressed in terms of the amplitudes to the even and odd scattering states, $\mathcal{A}_{oE\gets g}=\langle\psi_{oE}|\Psi_I(t)\rangle$ and $\mathcal{A}_{eE\gets g}=\langle\psi_{eE}|\Psi_I(t)\rangle$,
\begin{equation}
     \mathcal{A}_{L/R\, E\gets g}=\,-ie^{+i\eta_e}\,\mathcal{A}_{e E\gets g}\,\pm\, ie^{+i\eta_o}\,\mathcal{A}_{oE\gets g}.
\end{equation}
The even (odd) amplitudes result from the exchange, from the ground state, of an even (odd) number of photons. To lowest order, therefore, the left and right ionization amplitudes have the simplified form
\begin{equation}
     \mathcal{A}_{L/R\, E\gets g}\simeq\,-ie^{+i\eta_e}\,\mathcal{A}_{e E\gets g}^{(2)}\,\pm\, ie^{+i\eta_o}\,\mathcal{A}_{oE\gets g}^{(1)}.
\end{equation}
To the same order, the total ionization probability is
\begin{equation}
\mathcal{P}_{tot}(E)\simeq|\mathcal{A}^{(2)}_e(E)|^2+|\mathcal{A}^{(1)}_o(E)|^2,\label{eq:ionization_prob}
\end{equation}
whereas the left-right asymmetry is
 \begin{equation}
\Sigma(E) \simeq -4\,\Re e\,[\mathcal{A}_o^{(1)*}(E)\,\mathcal{A}_e^{(2)}(E)\,e^{i(\eta_e-\eta_o)}],
 \label{eq:asym}
\end{equation}
where $\Re e\, z$ is the real part of $z$. Equations~\eqref{eq:ionization_prob} and \eqref{eq:asym} show that the interference between one- and two-photon ionization pathways manifests itself in the photoemission asymmetry only.

To compute $\Sigma(E)$ perturbatively in the AM, it is necessary to evaluate the two-photon ionization amplitude \eqref{eq:two-photon-matrix}.
The two-photon matrix element in velocity gauge,
\begin{eqnarray}
\mathcal{M}_{eE\gets g}^{(2)}(\omega)&=&\mathcal{P}\int\frac{\langle \psi_{eE}|p|\psi_{oE'}\rangle\langle \psi_{oE'}|p|\psi_{g}\rangle}{E_g+\omega-E'}dE'-\nonumber\\ \label{eq:two-photon-matrix-2}
&-&i\pi{\langle \psi_{eE}|p|\psi_{o\,E_g+\omega}\rangle\langle \psi_{o\,E_g+\omega}|p|\psi_{g}\rangle},
\label{eq:two-photon-matrix-2}
\end{eqnarray}
has a principal-part $\mathcal{P}$ integral that we compute numerically from the dipole matrix elements, which are analytically known. The first term on the right hand side of \eqref{eq:two-photon-matrix-2}, which accounts for the contribution of the virtual excitations in the two-photon transition, is purely real, whereas the second term, which accounts for the sequential part of the excitation, is purely imaginary. Therefore, the phase $\varphi^{(2)}\approx\arg\left[\mathcal{M}_{eE\leftarrow g}^{(2)}(\omega_{\scriptscriptstyle\textsc{IR}})\right]$ depends on the ratio between the sequential and virtual components of the transition. Notice how in the soft-photon approximation~\cite{Galan16}, which is sometimes used to model multiphoton processes, the ionization continua are not coupled \emph{off shell}, and hence only the principal part contributes. In general, however, the continuum-continuum coupling does have off-diagonal terms. In the present case, the most relevant \emph{off-shell} contribution is due to the intermediate resonance. Finally, the two-photon ionization amplitude is computed using \eqref{eq:two-photon} and the left-right asymmetry is subsequently obtained using \eqref{eq:asym}.

\section{Results and discussion} \label{sec:Results}
In this section we discuss the results of the 1-2 RABBITT pump-probe ionization scheme applied to the ground state of both the analytical~\eqref{eq:pot_analy} and the numerical~\eqref{eq:pot_num} model. In both cases the two repulsive barriers are set at $10$~a.u. from the origin ($a=10$~a.u.). The potential strengths are chosen as $V_+=2.5$~a.u. and $V_-=0.5$~a.u., in the analytical model [see \eqref{eq:pot_analy}], whereas in the numerical model they are set to $V^N_+=3.1105$~a.u. and $V^N_-=0.5$~a.u. [see \eqref{eq:pot_num}], with a width of the Lorentzian potentials of $\Delta=0.05$~a.u. With these choices, the two models have similar ground state energies, $E_g=-0.125$~a.u. (AM) and $E_g=-0.1112$~a.u. (NM), as well as comparable positions and widths of the resonant states.

\subsection{Scattering phase and dipole matrix}
As discussed in the previous section, the scattering states of the AM are analytically known. To illustrate their essential features and compare them with the properties in the NM, we plot in Fig.~\ref{fig:phase-shifts} the total resonant phaseshift for AM and NM, from which we have subtracted the background phase $\eta_{\mathrm{bg}}=-2ka$, to highlight the resonant features: $\eta_{\mathrm{tot}}=\eta_e+\eta_o-\eta_{\mathrm{bg}}$. The resonant phase shift encodes the information needed to predict how the motion of a wavepacket in the potential differs from a reference free wavepacket. The derivative of $\eta_{\mathrm{tot}}(E)$ with respect to the energy is proportional to the averaged group delay experienced by a wave packet as it is transmitted or reflected by the potential. The steep $\pi$ jumps, in particular, indicate a long confinement, and hence they reveal the presence of resonant states. Due to the difference in the definition of the potential, the phase shifts $\eta_{e/o}(k)$ in the NM are not expected to exactly coincide with those in the AM. Nevertheless, as the figure shows, they do compare well with each other.
\begin{figure}[t]
\includegraphics[width=\columnwidth]{fig-3-phs-sft.eps}
\caption{Scattering phase shift for AM (blue solid line) and NM (red dotted line). Also shown the position of the two lowest odd resonant state.}
\label{fig:phase-shifts}
\end{figure}

The prominent role of resonances in the ionization of the ground state, in the present model, and the good agreement between the AM and NM, is apparent in the dipole transition matrix element $\mu_{E,g}$ [compare with Eq.~\eqref{eq:one-photon-mx}], which is shown in Fig.~\ref{fig:bc-dipole-matrix}.
\begin{figure}[b]
\includegraphics[width=1.\columnwidth]{fig-4-grd-sct-dtm}
\caption{\label{fig:bc-dipole-matrix} Square module of the bound-free dipole matrix element, as a function of the photoelectron energy calculated in AM (blue solid line) and NM (red dotted line).}
\end{figure}
Since the position of the resonances in the two models do not exactly coincide, the NM calculations in the three panels of Fig.~\ref{fig:bc-dipole-matrix} have been shifted to lower energies by $0.0015$~a.u., $0.0029$~a.u., and $0.0039$~a.u., respectively, to better highlight the similarity of the resonant peaks' shape and amplitude in the two calculations. Notice that the scale differ in each panel: resonance widths increase with energy.

\begin{figure}[t]
\includegraphics[width=1.\columnwidth]{fig-5-cnt-cnt-dtm.eps}
\caption{\label{fig:cc-dipole-matrix} Free-free transition dipole moments as a function of the even-state and odd-state energies for (a) AM and (b) NM.}
\end{figure}
Regarding the free-free transition dipole moments~\eqref{eq:ccdme}, we have verified that those computed in the NM are independent of the specific choice of $x_M$ in~\eqref{eq:boundaryCorrection}, as long as $x_M$ lies outside the effective range of the numerical potential barrier. The free-free dipole matrix elements for the two models are shown in Fig.~\ref{fig:cc-dipole-matrix} and found to be in excellent agreement with each other. In the figure, we observe the clear signature of the odd and even resonances, as well as the characteristic $(E-E')^{-1}$ singularity.

\subsection{Ionization amplitudes}
Let us consider a pump-probe scheme with electric field
\begin{equation}
E(t)=E_{\scriptscriptstyle\textsc{XUV}}(t)+E_{\scriptscriptstyle\textsc{IR}}(t,\tau),
\end{equation}
where the pump, which comprises two consecutive XUV harmonics, overlaps with a delayed IR probe pulse,
\begin{eqnarray}
E_{\scriptscriptstyle\textsc{XUV}}(t)&=&E_{\scriptscriptstyle\textsc{XUV}}F_{\scriptscriptstyle\textsc{XUV}}(t)\left[\sin(\omega_{2n}t)+\sin(\omega_{2n+1}t)\right],\nonumber\\
E_{\scriptscriptstyle\textsc{IR}}(t,\tau)&=&E_{\scriptscriptstyle\textsc{IR}}F_{\scriptscriptstyle\textsc{IR}}(t-\tau)\sin[\omega_{\scriptscriptstyle\textsc{IR}}(t-\tau)].
\end{eqnarray}
The two harmonics have equal maximum field amplitude $E_{\scriptscriptstyle\textsc{XUV}}=0.004$~au, a Gaussian envelope $F_{\scriptscriptstyle\textsc{XUV}}(t)$ with ${\rm FWHM}\approx 12$~fs, and $\Delta \omega=\omega_{2n+1}-\omega_{2n}=\omega_{\scriptscriptstyle\textsc{IR}}$.
The second harmonic $\omega_{2n+1}$ is tuned to the second odd resonance $E_1$ (see Fig.~\ref{fig:phase-shifts}), in both the AM and NM calculations.
The probe IR has a frequency $\omega_{\scriptscriptstyle\textsc{IR}}=0.051$~au, with ${\rm FWHM}\approx 57$~fs, and is delayed
with respect to the maximum of the harmonics' envelope by a time-delay $\tau$. The XUV and IR pulses are shown in Fig.~\ref{fig:pulse}.
 \begin{figure}[t]
\includegraphics[width=1.\columnwidth]{fig-7-pulse-num.eps}
\caption{\label{fig:pulse} XUV and IR pulses used in the study. Positive time delay $\tau$ indicates that the XUV harmonics arrive before the peak of the IR envelope.}
\end{figure}
\begin{figure}[b]
\includegraphics[width=1.\columnwidth]{fig-6-tot-specs}
\caption{\label{fig:tot-specs} Photoelectron distribution when using the XUV pump pulses only: AM (solid line), NM (dashed line). The ionization probability near the lower harmonic (grey region) has been multiplied by a factor of $1000$ to make it visible on the scale of the upper harmonic.}
\end{figure}

According to Eq.~\ref{eq:one-photon}, the one-photon ionization amplitude from the ground state, due to the XUV harmonics, is
 $ \mathcal{A}_{oE_f\leftarrow g}^{(1)} = -i\tilde{A}(\omega_{fg})\,\mu_{E_f,g}$, where $E_f$ is the final electron energy.
The ionization probability, $\mathcal{P}(E_f)=|\mathcal{A}_{oE_f\leftarrow g}^{(1)}|^2$, computed for the AM and the NM in the absence of the IR field are in excellent agreement, as shown in Fig.~\ref{fig:tot-specs}.
The ionization probability at the upper harmonic (resonance), induced by the higher harmonic $\omega_{2n+1}$, is orders of magnitude larger than that at the lower off-resonance harmonic $\omega_{2n}$, indicating that the one-photon background photoemission is negligible. This circumstance is to be expected since the ground state wavefunction, which is confined inside the potential barrier, has nearly zero overlap with the off-resonance free functions, which barely penetrate inside the barrier. This is also the reason why the resonances in our calculations have the symmetric Lorentzian shape characteristic of systems with negligible direct-ionization cross sections.
If the product $a\,V_-$ were much smaller, the ground state would eventually extend outside the confinement barrier, and the resonances  exhibit an asymmetric lineshape analogous to the celebrated Fano profile~\cite{Fano1961}.
\begin{figure}[b]
\includegraphics[width=1.\columnwidth]{fig-8a-asym-res}
\caption{\label{fig:8-res} Left-right asymmetry $\Sigma(E)$ as a function of the time-delay $\tau$ and the photoelectron energy, computed for AM (PT) (a) and NM (TDSE) (b), near the resonance $E_1$ where the dashed line indicates the center of the upper harmonic for each model.}
\end{figure}

\subsection{Two-Photon Ionization amplitude and asymmetry}
The left-right asymmetry \eqref{eq:asym} depends on the phases of the one- and two-photon ionization amplitudes
\begin{eqnarray}
\arg[\mathcal{A}_o^{(1)}(E)]&=&\varphi_{\scriptscriptstyle\textsc{XUV}},\\
\arg[\mathcal{A}_e^{(2)}(E)]&=&\varphi_{\scriptscriptstyle\textsc{XUV}}+\omega_{\scriptscriptstyle\textsc{IR}}\tau+\varphi^{(2)},
\end{eqnarray}
where $\varphi_{\scriptscriptstyle\textsc{XUV}}$ is the XUV phase, which is constant across the two harmonics, and $\varphi^{(2)}$ is the two-photon phase. The particular value of the relative phase between the two harmonics, which may depend on the process used to generate them, is irrelevant. For an IR pulse covering many cycles,
as it is the case in our study, $\varphi^{(2)}\approx\arg\left[\mathcal{M}_{eE\leftarrow g}^{(2)}(\omega_{\scriptscriptstyle\textsc{IR}})\right]$.
Therefore, the asymmetry~\eqref{eq:asym} beats at the same frequency as the IR field, as a function of the time delay $\tau$,
\begin{equation}
\Sigma(E) =\Sigma_0(E)\cos(-\omega_{\scriptscriptstyle\textsc{IR}}\tau+\varphi^{(2)}+\eta_e-\eta_o),
\label{eq:asym-dev}
\end{equation}
where $\Sigma_0(E) = -4|\mathcal{A}^{(1)}_{oE\gets g}\mathcal{A}_{eE\gets g}^{(2)}| $.
\begin{figure}[t]
\includegraphics[width=1.\columnwidth]{fig-8b-asym-sb}
\caption{\label{fig:8-sid} Same as in Fig.~\ref{fig:8-res} but at the lower harmonic and the dashed line now indicates the center of the lower harmonic for each model.}
\end{figure}
The NM photoionization asymmetry is readily obtained by projecting the final wavefunction $\Psi(t)$ from the TDSE simulation onto the left $\psi_{LE}$ and right $\psi_{RE}$ scattering states and using \eqref{eq:asymmetry}.
The left-right asymmetry, recorded as a function of the time delay between the XUV-pump and IR-probe, is plotted in Figs. \ref{fig:8-res} and \ref{fig:8-sid}, for both AM (PT) and NM (TDSE).
 The figures demonstrate the predicted $\omega_{\scriptscriptstyle\textsc{IR}}$ beating for each harmonic. In Fig.~\ref{fig:8-res} we present the asymmetry beating at the upper harmonic, near the resonance, whose phase experiences a dramatic excursion across the resonance. This result is easily understood from \eqref{eq:asym-dev}: whereas the phase $\varphi^{(2)}+\eta_e$ varies smoothly across the resonance, $\eta_o$ undergoes the same rapid excursion visible already in Fig.~\ref{fig:phase-shifts}. From Fig.~\ref{fig:8-res}, therefore, it is possible to extract the pure one-photon scattering phase associated with that resonant state.
 In contrast to the asymmetry beating at the resonant harmonic, the phaseshift at the lower harmonic in Fig.~\ref{fig:8-sid} shows a much less pronounced modulation, similar to the one that is accessible with traditional RABBITT measurements~\cite{Galan16,Kotur2016,Gruson2016}.
 The main phase variation is now due to $\varphi^{(2)}$, which encodes the resonance phase through the first XUV photon absorption, distorted by the free-free transition from the upper to the lower harmonic.
 Even in this case the PT and TDSE are in excellent agreement.
\begin{figure}
\includegraphics[width=0.8\columnwidth]{fig-9a-ext-comp-plane}
\caption{\label{fig:9} Averaged Fourier transform of the left-right asymmetry, with respect to the time delay, $\bar{\Sigma}(E)$ in the complex plane at the resonance (large circles, red) and at the lower harmonic (small circles, blue), for AM (solid lines) and NM (dashed lines). The arrows indicate the evolution of the trajectory with increasing energy.}
\end{figure}
To extract the phaseshift from the asymmetry beating at the upper harmonic we computed the Fourier transform
\begin{equation}
 \tilde{\Sigma}(E;\omega)=\frac{1}{\sqrt{2\pi}}\int^{\infty}_{-\infty}\xi(\tau)\,\,\Sigma(E;\tau)\,e^{-i\omega\tau}\,d\tau,
\end{equation}
where $\xi(\tau)$ is a smooth window function that allows us to limit the overlap between the FT of the background, centered at $\omega = 0$, with that of the beating signal, which is centered at $\omega=\omega_{\scriptscriptstyle\textsc{IR}}$,
\begin{eqnarray}
\xi(\tau)=\frac{1}{2}\left[1-{\rm erf}(2\tau+T_\xi)\,{\rm erf}(2\tau-T_\xi)\right],
\end{eqnarray}
${\rm erf(x)}=2/\sqrt{\pi}\int_0^x\exp(-t^2)dt$, and $T_{\xi}=9T_{\scriptscriptstyle\textsc{IR}}$, with $T_{\scriptscriptstyle\textsc{IR}}=2\pi/\omega_{\scriptscriptstyle\textsc{IR}}$ being the IR period.
Finally, $\tilde{\Sigma}$ is integrated across a spectral width $\Delta\omega\approx0.02$~a.u. around $\omega_{\scriptscriptstyle\textsc{IR}}$,
\begin{eqnarray}
\bar{\Sigma}(E)=\int_{\omega_{\scriptscriptstyle\textsc{IR}}-\Delta\omega/2}^{\omega_{\scriptscriptstyle\textsc{IR}}+\Delta\omega/2}\tilde{\Sigma}(E;\omega)\,d\omega.
\end{eqnarray}
The resulting complex amplitude, $\bar{\Sigma}(E)$, are plotted for both models in Fig.~\ref{fig:9}.
The complex amplitudes exhibit the characteristic jump of $\pi$ as the energy crosses the resonance while at the sideband the jump of $\pi$ is modulated by the emission of an IR-photon from a pulse of finite length.
\begin{figure}
\includegraphics[width=\columnwidth]{fig-10-ext-phs}
\caption{\label{fig:ion-asymm} Retrieved scattering phase $\eta_o^{ret} $obtained from the computed phase of $\bar{\Sigma}(E)$ (red solid line) and scattering phaseshift  $\eta_o$ (blue dashed line) near the resonance plotted for AM (PT) in (a) and NM (TDSE) in (b).  The asymmetry amplitude $|\bar{\Sigma}(E)|$ obtained from both models is also shown (green dashed line).}
\end{figure}
The resulting phases and amplitudes are plotted separately in Fig.~\ref{fig:ion-asymm}. We also compare the scattering phase at resonance retrieved from $\bar{\Sigma}(E)$ with the scattering phases $\eta_o$. To highlight the close similarity between the scattering phase and the phase retrieved from the pump-probe spectrum we have shifted the latter by a constant $\eta_c$. The striking agreement between the two phases demonstrates that the use of consecutive harmonics of the probe pulse holds the key to direct measurement of Wigner delays~\cite{Wigner1955,Pazourek2015}.
To extract the resonance position and width, the phase shift $\eta_{\mathrm{tot}}(E)$ is fitted around the resonance to the argument of the Breit-Wigner amplitude and a quadratic background $\delta_{\mathrm{bg}}(E)$ of the form:
\begin{equation}
\eta_{\mathrm{tot}}(E)=\acot[2(E-E_r)/\Gamma]+ \delta_{\mathrm{bg}}(E).
\end{equation}
The resonance width, $\Gamma$, and position, $E_r$, obtained from the scattering phaseshift are listed in Tab. \ref{tab:1} along side $\Gamma^{ret}$ and $E^{ret}_r$ fitted from $\eta^{ret}_o$.
\begin{table}[hbtp!]
\caption{\label{tab:1} Comparison of resonance position and width obtained directly from $\eta_o$ or from the retrieved  phaseshift  $\eta^{ret}_o$. All values are given in atomic units.}
\begin{ruledtabular}
\begin{tabular}{ccccc}
&$E_r$&$\Gamma\times10^3$&$E^{ret}_r$&$\Gamma^{\rm ret}\times10^3$\\
\hline
AM & 0.1913 & 1.184 & 0.1913 & 1.181 \\
NM & 0.1942 & 1.208 & 0.1942 & 1.183 \\
\end{tabular}
\end{ruledtabular}
\end{table}

\section{Conclusion} \label{sec:CONCLUSION}
We have presented 1-2 RABBITT, a new perturbative attosecond interferometric spectroscopy  in which two consecutive isolated harmonics of the probe pulse are used, instead of full comb of odd harmonics, as it is the case in traditional RABBITT. Using a 1D potential that supports bound as well as autoionizing states, we have shown how the new technique allows us to retrieve full phase information of the photoelectron-emission amplitude without the measurement-induced distortion inherent to traditional RABBITT. This feature results from the direct interference between one- and two-photon transition pathways. Furthermore, 1-2 RABBITT is also able to quantify the measurement-induced distortion due to the exchange of a probe photon. These findings have general validity and can be applied to study autoionizing states in real systems, such as rare-gas atoms. Specific generation schemes in which the APT comprises isolated consecutive harmonics of the probe, such as when the probe is the second harmonic of the fundamental frequency, or when the APT is generated with bicircular light pulses, as well as their application to real atomic systems, is beyond the scope of the present work and will be the subject of future investigations.

\section*{ACKNWLEDGMENTS}
This work was supported by the United States National Science Foundation under NSF Grant N$^\circ$~PHY-1607588, and by the UCF in-house OR grant Acc. N$^\circ$~24089045.


\section{Benchmarks}\label{sec:benchmarking}

\begin{table*}[!tbh]
  \centering

  \begin{tabular}{ l  l l  l*{16}{c}  l }
  \toprule   \\[-1em]
  \rotatebox{90}{\scriptsize Train} & \rotatebox{90}{\scriptsize Test} & \rotatebox{90}{\scriptsize road} & \rotatebox{90}{\scriptsize sidewalk} & \rotatebox{90}{\scriptsize person}  & \rotatebox{90}{\scriptsize motorcycle} & \rotatebox{90}{\scriptsize bicycle} & \rotatebox{90}{\scriptsize car} & \rotatebox{90}{\scriptsize truck} & \rotatebox{90}{\scriptsize bus} & \rotatebox{90}{\scriptsize wall} & \rotatebox{90}{\scriptsize fence} & \rotatebox{90}{\scriptsize traffic sign} & \rotatebox{90}{\scriptsize traffic light} & \rotatebox{90}{\scriptsize pole} & \rotatebox{90}{\scriptsize building} & \rotatebox{90}{\scriptsize vegetation} &  \rotatebox{90}{\scriptsize sky} & \rotatebox{90}{\parbox{1cm}{\scriptsize mIoU of\\[-0.2em] common\\[-0.2em] labels}} \\
  \midrule \\[-2ex]
  CS & DS & 72 & 22 & 30  & 47& 10&  58 & 30 & 19 & 17 & 13 &  19 & 8 & 23& 32& 76& 68& 34   \\
  DS & CS & 81 & 26 & 74  &  34 & 55 & 85 & 16 & 17 & 21 & 24 & 25 & 21 & 47& 77& 90& 88 & 49   \\

  \midrule \\[-2ex]
  BD & ID & 83 & 0 &  38 & 44 & 2 & 52 & 21 & 13 & 0 & 0 & 0 & 0 & 36 & 42 & 83 & 94 & 32   \\
  ID & BD & 84 & 16 &  57 & 34 & 44 & 77 & 14 & 24 & 10 & 33 & 18 & 13 & 41 & 68 & 82 & 87 & 44  \\
  \midrule \\[-2ex]
  CS & CS & 98 & 84 & 81 & 60 & 76 & 94 & 56 & 78 & 49 & 58 & 77 & 67 & 62& 92& 92& 94& 76 \\

  MV & MV & 85 & 58 & 73 & 55 & 61 & 90 & 61 & 65 & 45 & 58 & 72 & 67 & 50 & 86 & 90 & 98 & 70 \\
  ID & ID & 92 & 68 & 73 & 80 & 42 & 89 & 79 & 78 & 64 & 45 & 60 & 38 & 58 & 75 & 90 & 97 & 70 \\
  BD & BD & 95 &  62 & 61  & 32  & 22  & 90 & 52 & 57 & 25 & 45 & 52 & 58 & 49 & 85 & 87 & 97 & 60 \\

  \bottomrule \\[-1.5ex]

  \end{tabular}
  \caption{\label{tab:dom-desc} The domain discrepancy between Cityscapes (CS) \cite{CS}, Mapillary Vistas (MV) \cite{Mappilary}, Berkeley Deepdrive (BD) \cite{BDD100K} Dataset and IDD~ (ID) using the DRN-D-38 Model \cite{drn}. Performance for only the common labels between the four datasets are used. First two rows compares the accuracy of a model trained on one of  IDD~ or Cityscapes and tested on the other dataset. As can be seen, IDD~ trained model can predict CS and BD labels, better than predictions of trained models of the corresponding datasets on IDD. The bottom four rows gives the performance of models in each of the datasets. IDD~dataset is harder than CS dataset and similar in hardness to MV on these 16 labels. BD is harder because i.) it has night scenes ii.) the images are take from a dash board cam, hence has reflections from inside the car as well as distortions like rain drops on the mirror. }
  \vspace{-1.5em}
  \end{table*}

\subsection{Control Experiments}
In Table \ref{tab:control-expts}, we provide the results of some control experiments which provide upper bounds for IoU scores for models giving predictions at a given factor of the input resolution. We first  downsample  the ground truth by a given factor and then upsample it to the original image size for evaluation of average IoU at original scale. We provide the mean IoU scores of different levels of the hierarchy,  confirming  that low-resolution  processing contributes significantly to overall degradation of segmentation results.

\subsection{Domain Discrepancy}
Domain discrepancy studies the quantitative shift in data distributions between datasets. To understand it, we train a DRN-D-$38$ (Dilated Residual Networks \cite{drn}) model in Cityscapes \cite{CS}, Mapillary \cite{Mappilary}, BDD100K \cite{BDD100K} and our dataset. We compare the IoU scores in a set of 16 common labels between the four datasets in Table \ref{tab:dom-desc}. As seen from the last four rows, our dataset is harder than Cityscapes while having a similar level of hardness compared to the Mapillary Vistas Dataset. We also report IoU scores of predictions given by models trained on one dataset and tested on the other. A pretrained model trained on our dataset performs better when tested on Cityscapes and BDD100K, as compared to the converse experiment.

\subsection{Semantic Segmentation Benchmark}

\begin{figure*}
  \centering
  \includegraphics[width=\textwidth]{ious.pdf}
  \vspace{-1.8em}
  \caption{\label{fig:ious} The IoUs for every class for the DRN D 38 model trained on IDD with mIoU of 66.5\%.}
  \vspace{-1.3em}
  \end{figure*}
The semantic segmentation benchmark on our dataset quantifies the mean Intersection over Union (mIoU) scores at the four levels of the hierarchy. There are some labels in level 4 like traffic light, parking or animal for which the number of labeled pixels are very few. Hence, this serves as an excellent benchmark for transfer learning or domain adaptation problems. We also have level 1 and level 2 mIoU scores, which are a useful benchmark for real-time models, since they might not be able give good results on the fine grained classification task at level 3 and 4. The level 1 benchmark still has classes for \MC{most of} the essential labels for autonomous navigation.

We benchmarked our dataset using the DRN-D-38 (\cite{drn}) and the ERFNet (real-time model \cite{ERFNet}) model. We also conducted a challenge and evaluated submissions which use some of the state-of-the-art models. The results are shown in Table \ref{tab:sem-bench}.

\begin{table}
  \centering
  \begin{tabular}{l c c c }
    \toprule  \\[-1em]
    \multirow{2}{*}{Method} &
      \multicolumn{3}{c}{\parbox{2.5cm}{\centering \% mIoU at Levels} }  \\
      \cmidrule{2-4}
    & L1 &  L2 & L3  \\
    \midrule \\[-1.5ex]
    ERFNet & - & - & 55.4  \\[0.8em]

    DRN-D-38 & 85.9 & 72.6  & 66.6  \\[0.8em]

     *DeeplabV3+ \cite{DeeplabV3P} & 89.8 & 78.0  & 74.0 \\[0.8em]
     \parbox{4.8cm}{*PSPNet \cite{PSPNet}} & 89.9 & 78.0  & 74.1 \\[1.2em]

    \parbox{4.8cm}{*Wider Resnet-38, DeeplabV3 Decoder, Inplace ABN \cite{inplace}, Ensemble of 4} & 89.7 & 77.9  & 74.3 \\
    \bottomrule \\[-1.5ex]
  \end{tabular}
  \caption{\label{tab:sem-bench} The mIoU scores of models at 3 level of the hierarchy. The performance numbers of * models are obtained from the submissions of a AutoNUE challenge \cite{autonue} conducted based on the dataset.}
  \vspace{-1.3em}
\end{table}

\subsection{Class IoUs and Confusion Matrix}
The IoUs for every class can be found in Figure \ref{fig:ious}. We observe that the IoUs are lower than 25\% for bicycle, traffic light, vehicle fall-back and fence labels. The low scores for bicycle and traffic light can be explained by the low pixel counts.
We also plot the confusion matrix between labels in Figure \ref{fig:confusion}. Note that there is significant confusion between:
\begin{tight_itemize}
\item motorcycle and bicycle.
\item billboard and traffic sign.
\item obs-str-bar-fallback, vegetation and traffic light.
\item building and billboard.
\item vegetation and wall, pole, fence.
\item drivable, non-drivable, vegetation.
\end{tight_itemize}

\begin{figure}
    \centering
    \includegraphics[width=0.48\textwidth]{cm.jpg}

    \caption{\label{fig:confusion} The confusion matrix of the trained model.}
   \vspace{-0.2em}

\end{figure}

We analyze some examples of predictions in Figure \ref{fig:qualitative}. As can be seen, model trained on our dataset gives prediction of much better quality in unstructured setting. It identifies the muddy areas which can be driven. New labels like autorickshaw, curb, billboard etc are getting identified.

\begin{table}
  \begin{tabular}{l l l l l l l}
  \toprule
    Method
    & AP &  AP@50\\
    \midrule \\[-1.5ex]
    \parbox{5cm}{*MaskRCNN \cite{MaskRCNN} with ResNet101} & 0.268 & 0.499\\[1em]
    \parbox{4.5cm}{*PANet \cite{PANet}} & 0.376 & 0.661 \\
    \bottomrule \\[-1.5ex]
  \end{tabular}
  \caption{\label{tab:inst-bench} The AP scores of models for instance labels. The performance numbers of * models are obtained from the submissions of the AutoNUE challenge \cite{autonue} conducted based on the dataset.}
  \vspace{-1.2em}
\end{table}

\subsection{Instance Segmentation Benchmark}
Similar to other datasets, we also specified an instance segmentation benchmark, where individual instances of the same label need to be segmented separately. The algorithms are required to predict a set of detections of traffic participants in the frame,  with a confidence score and a per-instance binary segmentation mask.  To assess  the
performance, the average precision on the region level for each class and average it across a range of overlap thresholds ranging  from $0.5$ to $0.95$ in  steps  of  $0.05$, similar to \cite{coco}.

The results of some best performing submissions from the challenge are given  in Table \ref{tab:inst-bench}.

\begin{figure*}
\centering
\includegraphics[width=\textwidth]{autos.jpg}
\includegraphics[width=\textwidth]{shadow.jpg}
\includegraphics[width=\textwidth]{no_rules.jpg}

\caption{\label{fig:qualitative} We give many qualitative example with: input image from validation set, predictions from Cityscape pretrained model, prediction from model trained on our training dataset and the ground truth in our dataset in the order of columns.}
\end{figure*}

\section{Challenges in Unstructured Environments}
\label{sec:challenges-india}
We collect data from Indian roads and analyze the shortcomings of models trained on existing datasets. As illustration, we describe some of the qualitative issues observed when using predicted outputs of a model that obtains 70\% mean IoU on the Cityscapes validation set.

\vspace{-0.3cm}
\paragraph{Ambiguous Road Boundaries.}
Road boundaries in Cityscapes are very well defined and usually flanked on both sides by barriers or sidewalks. However, this is not the case in our setting. Road sides can have muddy terrain, while also being drivable to some extent. Roads themselves can be covered by dirt or mud, making the boundaries very ambiguous. On the other hand, Cityscapes models often recognize flat areas beside the road which need not be safe for driving as road, as  seen in Figure \ref{fig:roadside}.
\begin{figure}[!!t]
\centering
\includegraphics[width=0.49\textwidth]{124_roadside.png}
\includegraphics[width=0.49\textwidth]{124_roadside_2.png}
\vspace{-1em}
\caption{\label{fig:roadside} Cityscape models do not distinguish between the road and possible unsafe drivable area on both sides of the road.}
\vspace{-1em}
\end{figure}

\vspace{-0.3cm}
\paragraph{Diversity of Vehicles and Pedestrians.}
Indian roads have a variety of unique vehicles like auto-rickshaws, which behave very differently than other vehicles like cars. \MC{Even for standard categories like cars, the appearance variations are higher due to greater wear and tear.} Further, the frequency and variety of trucks and buses are also high. Another distinction is the large number of motorbikes with multiple persons riding it. Pedestrians often cross the road at arbitrary locations, rather than crosswalks. Bikes and autorickshaws \MC{are also less likely to follow traffic discipline, thus, there are fewer correlation between traffic participants and road signage such as lanes or traffic lights.}

\begin{figure}
\centering
\includegraphics[width=0.49\textwidth]{106_buffello_cs.png}
\includegraphics[width=0.49\textwidth]{auto_veh_120.png}
\includegraphics[width=0.49\textwidth]{auto_veh_120_2.png}
\includegraphics[width=0.49\textwidth]{bike_gas.png}
\vspace{-1em}
\caption{\label{fig:div-vehicles} (Top) A herd of buffaloes on the road at dusk. (Bottom) Many motorbikes with multiple riders, not necessarily following traffic rules such as road lanes.}
\vspace{-2em}
\end{figure}

\vspace{-0.3cm}
\paragraph{Extensive Use of Information Boards.}
Information displays such as billboards appear extensively in our dataset. They can be significant for localization and mapping problems by indicating buildings or landmarks. Sometimes they also indicate special vehicles, such as advertisements attached to a driving school car, or a delivery vehicle.

\begin{figure}
\centering
\includegraphics[width=0.49\textwidth]{130_info_boards.png}
\includegraphics[width=0.49\textwidth]{130_driv_school.png}
\vspace{-1em}
\caption{\label{fig:billboards} (Left) An array of billboards indicating the shops. (Right) A vehicle with a billboard of a driving school.}
\vspace{-2em}
\end{figure}

\vspace{-0.3cm}
\paragraph{\MC{Diversity of Ambient Conditions.}}
\MC{Lighting variation in our dataset is high since we acquire images at various times of the day, including mid-day, dawn and dusk.}
Also, some of the images have heavy shadows, which are common during a long summer season. Our dataset also contains scenes with heavily clouded skies. \MC{The greater variation in particulate matter due to fog, dust or smog also leads to significant appearance variations.} Cityscapes pretrained models yield lower accuracies in these settings, as seen in Figure \ref{fig:shade1}.

\begin{figure}
\centering
\includegraphics[width=0.49\textwidth]{shade1.png}
\includegraphics[width=0.49\textwidth]{low_light.png}
\vspace{-1em}
\caption{\label{fig:shade1} Heavy shadows in the image (top) or low light conditions (bottom) can greatly degrade the quality of predictions using models trained on Cityscapes. }
\vspace{-1em}
\end{figure}

\section{Conclusion}
We present a novel dataset for studying problems of autonomous navigations in unstructured driving conditions. We identify several drawbacks of existing datasets, such as distinguishing safe or unsafe drivable areas beside the road, additional labels required for vehicles and a label hierarchy that reduces ambiguity. We analyze the label statistics and the class imbalance present in the dataset. We also examine the domain discrepancy properties with respect to other semantic segmentation datasets. \MC{In contrast to existing datasets on semantic segmentation, ours is acquired in India, which leads to greater diversity due to variations in appearance of traffic participants as well as background categories. Not only does this pose interesting challenges for the state-of-the-art in semantic segmentation, it is also the first effort in our knowledge to focus on problems related to autonomous driving in geographies outside North America or Europe with relatively less developed road infrastructure.}

\MC{In the future, we plan to extend the benchmark to computer vision problems beyond semantic segmentation. In particular, the unconstrained nature of the dataset provides a uniquely novel setting for higher-level reasoning problems such as scene understanding \cite{su1, su2, su3} and path planning \cite{pp}. Motions in the dataset are less constrained due to greater freedom in traffic participant behavior and less adherence to traffic rules. The possible absence of visual cues such as lanes that constrain traffic participant behavior  poses further challenges. Besides the presence of rare categories \cite{rare}, even common categories have diverse attributes or appearance variations. Besides, the ambient conditions differ greatly across weather, time of day and air quality. This also motivates interesting new problems for few shot learning \cite{os} and domain adaptation \cite{da}, which our future work will study in greater detail.}

{ \small
\textbf{Acknowledgements.}
The authors would like to thank Intel, specially the Intel India
team for the efforts in capture of the data and coordination. Authors would like to specially acknowledge the support and
helps from Bharat Kaul, Prabhavathy Adavikolanu and Silar Shaik
in making this possible.
We would also like to thank Governments of Telangana and
Karnataka for the permissions, encouragement and enabling
this effort.
}
\section{Dataset}\label{sec:spec}

\begin{table*}
  \centering

    \begin{tabular}{l l l l l l l l l r }
    \toprule  \\[-1em]
   Dataset  & Calibration & \parbox{1cm}{Nearby\\frames\\ / Video} & \parbox{1.5cm}{Distortion\\ /Night}  &  \parbox{1.5cm}{\#Images/ \\ \#Sequences} & \parbox{1.5cm}{\#Labels\\Train/Total} & \parbox{1.5cm}{Average\\ Resolution}   \\[1.2em]
      \midrule \\[-2ex]
   Cityscapes \cite{CS}  & \checkmark &\checkmark & & 5K / ~50 & 19/34 & 2048x1024    \\
   IDD    & \checkmark & \checkmark & &  10K / 180 &  30/34  & 1678x968  \\
   BDD100K \cite{BDD100K}   &  & \checkmark & \checkmark  & 10K / 10K & 19/30 & 1280x720   \\
   MVD \cite{Mappilary}  & &  &  & 25K / ~-~~ & 65/66 & $>$1920x1080   \\

      \bottomrule \\[-1.5ex]
    \end{tabular}
    \caption{\label{tab:ds-comp} Comparison of semantic segmentation datasets for autonomous navigation.}

    \vspace{-1.2em}
  \end{table*}

\subsection{Acquisition}

The data was collected from Bangalore and Hyderabad cities \MC{in India} and their outskirts. The locations have a mix of urban and rural areas, highway, single lane and double lane roads with a variety of traffic. The driving conditions in these localities are highly unstructured due to multiple reasons: (i) these cities are rapidly growing and have a lot of construction near the roads, (ii) road boundaries are not well defined, (iii) pedestrians and jaywalkers are aplenty in these road images, and (iv) high density of  motorbikes and trucks on the road. The variety of vehicle models are also very large. A total of 182 drive sequences were used for the preparation of the dataset.

\subsection{Frame Selection}
 We chose images from one of the forward facing cameras of a stereo pair, for  fine annotation. Images were sampled at varying rates from the video sequence, with denser sampling around crowded and special interest places like traffic junctions. These images were annotated very finely, by layered polygon masks similar to Cityscapes. Since the road conditions are highly unstructured, we need a wider variety of labels. We annotated a total of 10,004 frames.

\begin{figure*}
  \centering

  \includegraphics[width=\textwidth]{pixel_cs_compare.pdf}
  \vspace{-1em}
  \caption{\label{fig:pix_comp_cs} Comparison of the pixel count in our dataset with Cityscapes. The y axis is plotted in log-scale. Note that for most classes of vehicles, the number of pixels are 5-10 times more than Cityscapes. Moreover our dataset has newer labels like autorickshaw, billboard, drivable/nondrivable fallback which also have significant number of labeled pixels.}
  \vspace{-1.2em}
  \end{figure*}

\subsection{Label Hierarchy \& Annotation}
We used a total of 34 labels in the fine annotations. The labels were given definition by means of a textual description as well as example images. However, we found that it is difficult to completely avoid ambiguity between some labels. For example labels like parking, caravan or trailer cannot be precisely defined due to the diversity of the scenes and vehicles in the data collected. For resolving this issue, we designed a 4 level label hierarchy having 7 (level 1), 16 (level 2), 26 (level 3) and 30 (level 4) labels (see Figure \ref{fig:label-heirarchy}). Each level defines a category as the union of labels in the succeeding level, which are chosen such that they are ambiguous. Since we take unions of the most ambiguous labels while designing the hierarchy, the lower levels have lesser ambiguity. We have a set of new labels not available in Cityscapes \cite{CS} like auto rickshaw, billboards, animal, curb. We also have separate labels for road, drivable fall-back and non-drivable fall-back indicating safe, unsafe and non-drivable flat surfaces. We have added fall-back labels whenever appropriate so that highly ambiguous objects can be given labels.

For labeling the dataset, the annotation team was first asked to re-annotate images from the Cityscapes \cite{CS} dataset. The difference between the annotations were subsequently shown to the annotators. This process was done until the annotators were achieving greater than $95$\% accuracy with respect to the Cityscapes ground truth labels.

\subsection{Statistical Analysis and Dataset Splits}
The pixel statistics among the labels can be seen in Figure \ref{fig:label-heirarchy}. The labels in level 4 have high class imbalance. Labels like parking, animal, caravan or traffic light have much fewer pixels. The annotated dataset also has labels for trailer and rail track, which were combined with vehicle fallback and nondrivable fallback in level 4, since they have very few pixels that mostly fell within a few drive sequences. The lower levels are designed such that the imbalance is lesser.

Class imbalance at level 4, creates a problem while splitting the dataset in train, test and validation sets. Since the splitting is done at the level of drive sequences (that is, all images within a drive sequence are moved to the same split), we need to be careful that the few drive sequences that contains a label are rightly split. We roughly divide the dataset in to 70\% train, 10\% validation and 20\% test splits. The split was done by randomly assigning the drive sequences with the said distribution. We did the splitting multiple times to come up with a split where all the 30 labels in level 4 have approximately 70, 10 and 20 percentage of pixels in  the train, validation and test sets, respectively.

\subsection{Comparison with Other Datasets}
\begin{figure}
  \centering
  \includegraphics[width=0.45\textwidth]{instances_cs_compare.pdf}
  \caption{\label{fig:inst_tp} Comparison of traffic participants in our dataset with Cityscapes.}
  \vspace{-1.5em}
\end{figure}
Various datasets have been proposed for studying the semantic and instance segmentation problems like Pascal VOC \cite{pascal}, MS COCO \cite{coco}, SUN \cite{sun,sun2}.
The datasets proposed for semantic segmentation that focus on autonomous navigation are Cityscapes \cite{CS}, KITTI \cite{KITTI}, Camvid \cite{camvid,camvid2}, Leuven \cite{leuven} and the Daimler Dataset \cite{dailmerdataset}. ADE20K \cite{ADE20K} is a recent dataset which focuses on the general scene parsing problem. More recently Mapillary Vistas dataset \cite{Mappilary} (which focuses on street view imagery) and the Berkeley Deep Drive Dataset \cite{BDD100K} (for autonomous navigation) was released. A comparison of the metadata available in these datasets can be seen in Table \ref{tab:ds-comp}. As can be seen, our dataset is more similar to Cityscape in the sense that we collect data from calibrated cameras without any distortions. BDD100K uses dashboard cameras kept inside the car and hence often has internal reflections from the glass as well as rain distortions. Moreover a good fraction of the dataset consists of night images. Mapillary Vistas dataset consists of images taken using a variety of cameras (including smart phones) having varying perspectives of the road and the road side. They do not have video data or images of near by frames.

\begin{table}
  \centering
  \begin{tabular}{l c c c}

      \toprule   \\[-1em]
      \multirow{2}{*}{\parbox{2.5cm}{Method}} &
        \multicolumn{3}{c}{\parbox{4cm}{\centering \% mIoU at {Hierarchy Levels}}}  \\
        \cmidrule{2-4}
      & \parbox{1.2cm}{\centering 1} &  \parbox{1.2cm}{\centering 2} &  \parbox{1.2cm}{\centering 3} \\
      \midrule \\[-2ex]

      GT subsampled by 2 & 99 & 97 & 97\\

      GT subsampled by 4 & 98 & 96  & 95\\

      GT subsampled by 8 & 96 & 92  & 90\\

      GT subsampled by 16 & 92 & 87  & 84\\

      GT subsampled by 32 & 86 & 78  & 74\\

      GT subsampled by 64 & 77 & 66  & 61\\

      GT subsampled by 128 & 65 & 53  & 48\\

      \bottomrule \\[-1.5ex]
    \end{tabular}
    \caption{\label{tab:control-expts}  Control  experiments  to  estimate  upper  bounds  for  se-
    mantic segmentation results, assessed by Intersection-over-Union
    (IoU, in \%) scores for different levels of the hierarchy.}
    \vspace{-1em}
  \end{table}

We compare the label statistics of our dataset with Cityscapes (since it is more similar to our dataset as described above) in terms of pixel counts (Figure \ref{fig:pix_comp_cs}) and the number of instances of traffic participants (Figure \ref{fig:inst_tp}). We have more pixels of truck, bus, motorcycle, guard rail, bridge and rider (see Figure \ref{fig:pix_comp_cs}). The pixel counts for new labels (auto rickshaw, billboard, curb, drivable-fallback, nondrivable-fallback) are also high. In terms of instances of traffic participants, we have almost double the counts, with a distribution similar to Cityscapes (see Figure \ref{fig:inst_tp}).

\section{Introduction}
\label{sec:intro}

Autonomous navigation is rapidly maturing towards becoming a mainstream technology, with even consumer deployment by major automobile manufacturers. A significant contributor to this progress has been the availability of large-scale datasets for sensing and scene understanding. Yet, several challenges remain in enabling self-driving across diverse geographies. A key challenge is to achieve data scale and diversity large enough to ensure safety and reliability in extreme corner cases. Even more importantly, algorithms are largely untested in their ability to generalize to road conditions that are significantly more diverse and unstructured.

In this paper, we propose IDD, a dataset that takes the first steps towards addressing the above concerns. Our dataset shares several traits such as scale, annotation and tasks with similar ones in structured environments, namely KITTI \cite{KITTI} or Cityscapes \cite{CS}. But it also intends to significantly expand the scope of the autonomous navigation problem, along each of those dimensions.
Similar to Cityscapes \cite{CS}, we provide large-scale raw data with multiple cameras and sensors across cities and lighting conditions. But the scale of our data is larger, consisting of $10,004$ labeled images with fine instance-level boundaries. Next, while the annotation type is similar for our dataset and Cityscapes \cite{CS}, the number of object classes and within-class diversity of appearance are higher for IDD. Finally, while we also initially propose instance segmentation as the task of interest, our label diversity and novel hierarchy might allow novel machine learning techniques and computer vision algorithms.

\begin{figure*}[!th]
    \vspace{-4em}
    \centering
    \includegraphics[width=1.\textwidth]{pixels_label_heirarchy_legend.pdf}
    \vspace{-3em}
    \caption{\label{fig:label-heirarchy} Label distribution in our dataset. The following information is shown here: (i) pixel counts of individual labels on the y-axis (ii) four-level label hierarchy used by the dataset at the bottom, (iv) the color legend for the predicted and ground truth masks shown in the paper is used for the corresponding bars.  We define metrics at 4 levels of the hierarchy with 30 (level 4), 26 (level 3), 16 (level 2) and 7 (level 1) labels, respectively, giving different complexity levels for training models. }
    \vspace{-1em}
\end{figure*}

The singular defining aspect of our dataset is that it corresponds to driving in less structured environments. We argue that this is a better reflection of the needs for autonomous navigation in large portions of the world, including Asia, South America and Africa. Accordingly, we collect our data in India where road scenes differ markedly from those in Europe or North America. The variety of traffic participants in Indian roads is larger, including novel classes such as autorickshaws or animals. The within-class diversity is also higher, for example, since vehicles span a larger range of manufacturing years and ply with larger variation in wear. Even the distribution of classes that overlap with Cityscapes is significantly different, for instance, the proportion of motorcycles is far higher, as is that of multiple riders on two-wheeled vehicles. Background classes also display greater diversity, such as city scenes rich in novel classes such as billboards. Besides variations in weather and lighting, other ambient factors such as air quality and dust also span greater ranges in our dataset. Such greater complexity of road scenes necessitates a larger scale of data. Thus, we provide high-quality annotation at a scale significantly larger than available for other contemporary datasets such as KITTI or Cityscapes.
\footnote{Datasets such as Berkeley Deep Driving \cite{BDD100K} and Apolloscape \cite{apolloscape}  have recently been released with labels at a similar scale. However, they are contemporaneous with our work, which precludes a detailed comparison. In any case, we note that they are in structured environments, which makes our dataset clearly different.}

We provide a detailed analysis of the label distributions in the IDD dataset, while highlighting some of the above differences. We also showcase those differences through quantitative evaluation of state-of-the-art algorithms on our dataset. We consistently observe that semantic segmentation performances are far lower on IDD as compared to Cityscapes, using identical models and with larger-scale training data for IDD. Firstly, this highlights that conventional semantic segmentation datasets such as Cityscapes are getting saturated and the next set of challenges lie in more complex datasets like IDD. Secondly, this highlights the need for ever-larger training data as we expand the scope of the autonomous navigation problem to newer geographies.

Besides segmentation, the nature of our dataset also enables novel problems for vision and learning. This is already reflected in some of our annotation choices. For instance, while the notion of a drivable area in Europe is largely defined by classes such as roads or lanes that have distinct appearances, it is more ambiguous in our dataset and likely also informed by semantic cues such as presence of dynamic traffic participants. Thus, we include labels for safely drivable and non-drivable areas. Our label hierarchy is attuned to the autonomous navigation problem and we postulate that exploiting it might lead to semantic segmentation more suited to subsequent applications such as collision avoidance or path planning. We label classes that are rare but important for navigation (such as animals), or classes that exhibit large within-class variance (such as autorickshaw), which motivates problems such as few-shot learning.

The contrast of our dataset with structured ones also suggests interesting directions of future research. For instance, domain adaptation between Cityscapes and IDD is clearly a need given the large performance drops encountered in cross-dataset settings. Classes that are unique to our dataset also encourage consideration of domain adaptation with non-overlapping label spaces. Even higher-level reasoning problems such as behavior prediction pose new challenges in IDD, since traffic participants have lower adherence to traffic rules, motions can be sudden, complex obstructions might be present, drivable areas can be ambiguous and traffic lanes need not correspond to lane markings on the road. While not considered in this paper, we highlight that these novel problems do arise in unstructured environments such as ours.


\section{Introduction}
\label{introduction}
Electronic Limit Order Books (LOBs) \cite{parlour2008limit} are commonly used in modern financial markets. We can consider LOBs as the most granular financial data that records all limit orders\footnote{Limit orders are resting order to either buy or sell an instrument, i.e. those that do not match immediately. All resting limit orders form the LOB which represents the supply and demand of a given financial instrument at any moment in time.} of an exchange. Limit order data, like any other financial time-series data, is generally non-stationary and has a notoriously low signal-to-noise ratio.

Many models used historically for LOB data, such as the vector autoregressive model (VAR) \cite{zivot2006vector}, operate on handcrafted features extracted of the data to improve predictability. It is far from trivial to select informative features by hand. Recent research \cite{tsantekidis2017forecasting, sirignano2018universal, zhang2018deeplob} has shown that deep learning models can deliver good predictive performance using only raw LOBs data with the process of feature extraction being automated by convolutional layers. However, such deep learning techniques are (typically) non-probabilistic and thus provide only heuristic uncertainty measures associated with output variables. For any decision making process, especially for high-frequency trading, trust and risk are primary concerns.

Uncertainty may be quantified through Bayesian inference. Given the complexity of network models, such \emph{Bayesian Neural Networks} \cite{bishop1995neural} are often achieved by approximation such as variational inference \cite{bayes2008variational}. The work in \cite{gal2016dropout} proposed \emph{dropout variational inference}, also known as dropout sampling, as an approximation to BNNs. While variational dropout has been well established and deployed to many application domains, we have not found any literature applying it to LOBs.

\emph{Our contributions:} We apply variational dropout to a large-scale deep learning model designed in \cite{zhang2018deeplob} to predict short term movements in cash equity prices from LOBs. We showcase that dropout techniques, as stochastic regularizers, achieve better predictive performance and we design trading strategies that can deliver improved returns using derived uncertainty information. Further, we show how uncertainty relates to position sizing, avoiding bad trades and amplifying profits.

\section{Methods}
\label{methods}
We here briefly introduce the dropout variational inference \cite{gal2016uncertainty} and the network architecture, known as DeepLOB, in \cite{zhang2018deeplob}. In addition, we discuss different trading strategies and demonstrate how uncertainty information derived from Bayesian networks can be used to decide entry and exit points and further to upsize or downsize our positions.

\subsection{Dropout Variational Inference and Network Architecture}

For Bayesian networks that have more than one hidden layer, the true posterior distributions of model parameters can not usually be evaluated analytically. Instead, we can define a variational distribution, which has a simple structure to evaluate, to approximate the true posterior distribution. Following this variational interpretation, the work of \cite{gal2016uncertainty} showed that dropout can be considered as an approximating distribution to the posterior in a Bayesian network. In other words, the weights obtained from the optimisation of a network using dropout are the same as the variational parameters of a network that has the same structure. We can then model uncertainty in Bayesian networks by a technique called \textit{Monte-Carlo (MC) dropout}, which is essentially running a number of stochastic forward passes through our network and averaging the results,

The network architecture in \cite{zhang2018deeplob} is composed of three main building blocks: convolutional layers, an Inception Module and LSTMs. It is often difficult to extract representative features in financial applications as there are enormous observations and financial data is notoriously noisy. Hand-engineered features, such as technical indicators or PCA, often rely on tacit assumptions that may not well describe financial data and the process is non-trivial. The usage of CNNs and Inception Modules is then to automate the process of feature extraction, removing the above constraints. Feature maps obtained from convolutional layers are then passed into a LSTM layer instead of fully connected layers since these features still contain time dependencies. We can use LSTMs to capture temporal relationship among the resulting feature maps.

Dropout was not used in the above work but we apply it to obtain uncertainty estimates. However, in order to get well-calibrated uncertainty, we need to adapt the dropout probability and this process is often time consuming for the large model and dataset. As a result, we apply a dropout layer with rate 0.2 after the Inception Module, decided using grid-search methods. Note that the feature maps are still time-series, so we keep the dropout mask the same for all timesteps. Doing this would not break time dependencies but regularize resulting features. Further, we applied Concrete Dropout \cite{gal2017concrete} to convolutional layers in the Inception Module and found that all dropout rates are pushed towards zero due to the large number of observations we have available. Further, we note that the convolutional layers have considerably less parameters compared to the dense layers because of the parameter sharing structure, so they are less likely to be overfitted \cite{wu2015towards}.

\subsection{Different Trading Strategies}

The work of \cite{zhang2018deeplob} classified future price movements into three classes (up=+1, neutral=0, down=-1) and used a softmax layer to give a probability vector at each time stamp ${\mathbf{p}}_t = [p_{+1,t}, p_{0,t}, p_{-1,t}]^T$. The basic trading strategy designed in \cite{zhang2018deeplob}  to test the quality of the prediction was based on discrete outputs and did not use the softmax probability vector. We briefly describe this strategy below and refer to it as \textbf{Normal trading strategy}. We then show how returns can be amplified by incorporating information from the softmax layer (\textbf{Softmax trading strategy}) and finally show how positions can be sized by using uncertainty information from the Bayesian networks (\textbf{Bayesian trading strategy}).

\paragraph{Normal trading strategy}

The model delivers a softmax probability vector, ${\mathbf{p}}_t$, at each time. We enter the long position if $max({\mathbf{p}}_t) = p_{+1,t}$ or the short position if $max({\mathbf{p}}_t) = p_{-1,t}$. Only one position is allowed and we hold until the opposite movement is predicted (i.e. $max({\mathbf{p}}_{t+k}) = p_{-1,t+k}$ in case of a long position where $k$ is the holding period). We go long or short (buy or sell) the same number of shares ($\mu$) at each time and set $\mu$ to be 30\% of first ask or bid level depending on the position we enter.

\paragraph{Softmax trading strategy}

Same as above, but we use softmax probability vector to constrain our timing of entering positions. Instead of entering the trades when $max({\mathbf{p}}_t)$ is equal to $p_{+1,t}$ or $p_{-1,t}$, we place a simple rule to only enter the trade if $max({\mathbf{p}}_t)$ is greater than a threshold $\alpha$. We place no constraint to exit the trade.

\paragraph{Bayesian trading strategy}

While we use softmax outputs in the previous strategy, the network outputs could be biased and it is often argued \cite{gal2016uncertainty} that softmax probability vectors can not be interpreted as true model confidences. Instead, three methods, variation ratios, predictive entropy and mutual information are suggested in \cite{gal2016uncertainty} to summarise uncertainty in classification when using variational dropout. We tested on all three quantities and discover that predictive entropy, $\mathbb{H}$, gives most consistent uncertainty information and we use this information to size our positions. For a given input $\bm{x}$, we run $M(=100)$ stochasitc forward passes through our network and average the probabilities of each class to obtain the averaged probability vector $\bm{\bar{{\mathbf{p}}}}_t$. As before, we still enter trades if $max(\bar{{\mathbf{p}}}_t) > \alpha$, but we upsize our positions to $1.5 \times \mu$  if $\mathbb{H}<\beta_1$ or keep original size $\mu$ if $\beta_1<\mathbb{H}<\beta_2$ but downsize our positions to $0.5 \times \mu$ if $\mathbb{H}>\beta_2$. We exit the position if $\mathbb{H} < \beta_2$. In next section, we show  how these trading parameters affect our profits and risk level.

\section{Results}
\label{results}
The data we use are limit order book data for five stocks: Lloyds Bank, Barclays, Tesco, BT and Vodafone listed on London Stock Exchange for the entire 2017 year, totalling more than 134 million market quotes. We take the first 6 months as training data, the next 3 months as validation data and the last 3 months as test data. Our test data contains millions of observations that can thoroughly verify our model performance. We follow \cite{zhang2018deeplob} to prepare our input data.

We denote our model as BDLOB and compare it to two baseline models: the CNN model presented in \cite{tsantekidis2017forecasting} and DeepLOB5 in \cite{zhang2018deeplob}. We show each model's precision, recall, F1 score and average AUC score in Table~\ref{table:metrics} and we present confusion matrices and boxplots for daily accuracy in Figure~\ref{fig:cm_daily_accuracy}. BDLOB achieves the best results since dropout, as stochastic regularizers, improves generalization ability. Also, from Figure~\ref{fig:cm_daily_accuracy}, we can see that BDLOB has better ability to classify stationary class and lead to consistent better performance with narrow interquartile range (IQR) and less outliers.

\begin{table}[h]
  \caption{Evaluation metrics for different baseline models}
  \label{table:metrics}
  \centering
  \begin{tabular}{l|llll}
    \toprule
    Model     						& Precision     & Recall  &  F1 Score & AUC \\
    \midrule
    CNN   								& 0.52   & 0.52  	 & 0.48    & 0.672\\
    DeepLOB5   					& 0.60  	& 0.60    & 0.58  & 0.803\\
    BDLOB   & 0.60  	& 0.61    & 0.60  & 0.811\\
    \bottomrule
  \end{tabular}
\end{table}

\begin{figure}[h]
  \centering

  \includegraphics[width=5.5in, height=2.5in]{confusion_boxplot.pdf}
  \caption{\textbf{Top:} confusion matrices for CNN, DeepLOB5 and BDLOB; \textbf{Bottom:} boxplots of daily accuracy for the different models across 5 stocks.}
  \label{fig:cm_daily_accuracy}
\end{figure}

We follow the same assumptions, mid-price simulation and no transaction costs, as in \cite{zhang2018deeplob} and present how different trading strategies affect profitability. Our focus is not to design a stand-alone trading strategy but to showcase how uncertainty information can be beneficial to decide enter and exit points and also position sizing. However, due to the short holding periods, most high-frequency trading firms are engaged in market making. Being able to trade passively on one side of the trade effectively corresponds to a mid-mid trade, being able to trade both sides passively reduces the fee further. The above, simplified trading simulation still approximates the relative value of the here presented model when added to an existing market making strategy.

We hold no stocks overnight and close all positions by the end of the day. As each strategy leads to different trading volumes, we normalise each day's return by total volumes executed on that day. In the literature, the Sharpe Ratio \cite{sharpe1994sharpe} is often used to measure a strategy's risk. However, the Sharpe Ratio is based on variance of returns that treats both positive and negative returns to be ``risky''. In practice, few traders would consider large positive returns as risky so we use the Downward Deviation Ratio (DDR) \cite{moody2001learning} as a measure of risk. The DDR penalizes  negative returns and rewards large positive returns. The DDR is defined as,
\begin{equation}
\mathrm{DDR} = \frac{E(R_t)}{\mathrm{DD}_{T}} \quad \text{where} \quad \mathrm{DD}_{T} = \Big(\frac{1}{T}  \sum_{t=1}^T \min(R_t, 0)^2 \Big)^{\frac{1}{2}}
\end{equation}
Figure~\ref{fig:boxplot_return} presents the boxplots for daily profits normalized by total volumes on that day. Any strategy using uncertainty information leads to substantial improvements in returns and BDLOB delivers the most profits and maintain a similar risk level.  Figure~\ref{fig:cum_return} shows the cumulative profits across testing periods and we observe consistent and better performance for BDLOB. We set $\alpha = 0.7$ for all strategies that use softmax probability vectors and set $\beta_1=0.1, \beta_2=0.9$ for Bayesian trading strategy. The value of these hyperparameters are selected to balance profits and risk. We discuss these choices below.

Figure~\ref{fig:boxplot_hyper} shows how profit and risk change according to these hyperparameters. We can see that a higher $\alpha$ generally improves profits but $\alpha=0.8$ or $\alpha=0.9$ can lead to larger risk or lower profits -- this is due to few trades taking place so any wrong prediction has a stronger impact on performance. As a result, we select $\alpha=0.7$ as the best value and apply it to all Bayesian trading strategies in Figure~\ref{fig:boxplot_hyper}. We set $\beta_1 = 0.1$ and find profit or risk is not hugely affected by this value because our model is very confident in its predictions when a high $\alpha=0.7$ is used. The effect of $\beta_2$ is obvious as it not only penalizes unnecessary trades but also controls the exit time point. We can see that high profits can be obtained with a smaller $\beta_2$ but with a much higher risk level. A strict exit rule allows us to take even fewer trades amplifying profits, but the impact of wrong predictions is also more obvious and we can hold a bad trade for too long without leaving it, leading to large variance.
\begin{figure}[h]
  \centering
  \includegraphics[width=5.3in, height=2.9in]{boxplot_profit_final.pdf}
  \caption{Boxplots and DDR for normalized daily profits. Profits are in GBX(=GBP/100).}
  \label{fig:boxplot_return}
\end{figure}

\begin{figure}[h]
  \centering

  \includegraphics[width=5.5in, height=2.6in]{cumulative_profit_final.pdf}
  \caption{Cumulative profits for test periods for different models. Profits are in GBX(=GBP/100).}
  \label{fig:cum_return}
\end{figure}

\begin{figure}[h]
  \centering
  \includegraphics[width=5.5in, height=3.9in]{boxplot_profit_hyperparameter.pdf}
  \caption{Boxplots for normalized profits and DDR for different trading parameters. We denote a softmax strategy with rate $\alpha$ as BDLOB\_softmax\_$\alpha$ and Bayesian strategy as BDLOB\_bayesian$[ \beta_1, \beta_2]$. We set $\alpha=0.7$ for all Bayesian trading strategies.  Profits are in GBX(=GBP/100).}
  \label{fig:boxplot_hyper}
\end{figure}

\section{Conclusion}
\label{conclusion}
We show that Dropout inference improves performance for predicting stock price movements and demonstrate that uncertainty information can be utilised to design better trading strategies. A promising direction for future work is to explore recent Bayesian techniques such as \cite{khan2018fast} and \cite{badam} and consider trading strategies that monitor uncertainty information during holding periods.

\input{bayesian_deeplob.bbl}

\end{document}