We study the long time behavior of an advection-diffusion equation with a random shear flow which depends on a stationary Ornstein-Uhlenbeck (OU) process in parallel-plate channels enforcing the no-flux boundary conditions. We derive a closed form formula for the long time asymptotics of the arbitrary $N$-point correlator using the ground state eigenvalue perturbation approach proposed in \cite{bronski1997scalar}. In turn, appealing to the conclusion of the Hausdorff moment problem \cite{shohat1943problem}, we discover a diffusion equation with a random drift and deterministic enhanced diffusion possessing the exact same probability distribution function at long times. Such equations enjoy many ergodic properties which immediately translate to ergodicity results for the original problem. In particular, we establish that the first two Aris moments using a single realization of the random field can be used to explicitly construct all ensemble averaged moments. Also, the first two ensemble averaged moments explicitly predict any long time centered Aris moment. Our formulae quantitatively depict the dependence of the deterministic effective diffusion on the interaction between spatial structure of flow and random temporal fluctuation. This noteably contrasts the white noise case from the OU case. Further, this approximation provides many identities regarding the stationary OU process dependent time integral. We derive explicit formulae for the decaying passive scalar's long time limiting probability distribution function (PDF) for different types of initial conditions (e.g. deterministic and random). All results are verified by Monte-Carlo simulations.