The finite volume method (FVM) uses cell averages to represent the numerical solution. For example, $x_i$ indicates the center of the $i$-th grid cell, and $\textbf{U}_i^n$ denotes the numerical solution inside grid cell $i$ at time $t^n$. The fluxes at cell boundaries are indicated by $i \pm 1/2$ subscripts, e.g., using your paper's notation, $\Phi_{i+1/2}$ is the numerical flux at the right boundary of the $i$-th cell (at grid coordinate $x_{i+1/2}$).
Your grid cell indices range from $1$ to $N$. In addition to these, the authors define two "ghost" cells with indices $0$ and $N+1$. These are not part of your "real" mesh configuration. Instead, their definition allows you to compute the fluxes at the "real" boundaries of your grid (i.e., $\Phi_{1/2}$ and $\Phi_{N+1/2}$). In 1D, you can think of cell no $0$ as being immediately to the "left" of grid cell $1$. Similarly, cell $N+1$ lies immediately to the "right" of cell $N$, like in the illustration given here.
The solution approximation inside the ghost cells is computed using standard extrapolation formulas.