  A framework for integrating transfer matrices with particle-in-cell
simulations is developed for TeV staging of plasma wakefield accelerators.
Using nonlinear transfer matrices in terms up to ninth order in normalized
energy spread $\sqrt{\langle\delta\gamma^2\rangle}$ and deriving a compact
expression for the chromatic emittance growth in terms of the nonlinear matrix,
plasma wakefield accelerating stages simulated using the three-dimensional
particle-in-cell framework OSIRIS 4.0 were combined to model acceleration of an
electron beam from 10 GeV to 1 TeV in 85 plasma stages of meter scale-length
with long density ramps and connected by simple focusing lenses. In this
calculation, we find that for initial relative energy spreads below $10^{-3}$,
energy-spread growth below $10^{-5}$ of the energy gain per stage and
normalized emittance below mm-mrad, the chromatic emittance growth can be
minimal. The technique developed here may be useful for plasma collider design,
and potentially could be expanded to encompass non-linear wake structures and
include other degrees of freedom such as lepton spin.
