  Line intensity mapping (LIM) is a rapidly emerging technique for constraining
cosmology and galaxy formation using multi-frequency, low angular resolution
maps. Many LIM applications crucially rely on cross-correlations of two line
intensity maps, or of intensity maps with galaxy surveys or galaxy/CMB lensing.
We present a consistent halo model to predict all these cross-correlations and
enable joint analyses, in 3D redshift-space and for 2D projected maps. We
extend the conditional luminosity function formalism to the multi-line case, to
consistently account for correlated scatter between multiple galaxy line
luminosities. This allows us to model the scale-dependent decorrelation between
two line intensity maps, a key input for foreground rejection and for
approaches that estimate auto-spectra from cross-spectra. This also enables LIM
cross-correlations to reveal astrophysical properties of the interstellar
medium inacessible with LIM auto-spectra. We expose the different sources of
luminosity scatter or "line noise" in LIM, and clarify their effects on the
1-halo and galaxy shot noise terms. In particular, we show that the effective
number density of halos can in some cases exceed that of galaxies,
counterintuitively. Using observational and simulation input, we implement this
halo model for the H$\alpha$, [Oiii], Lyman-$\alpha$, CO and [Cii] lines. We
encourage observers and simulators to measure galaxy luminosity correlation
coefficients for pairs of lines whenever possible. Our code is publicly
available at https://github.com/EmmanuelSchaan/HaloGen/tree/LIM . In a
companion paper, we use this halo model formalism and code to highlight the
degeneracies between cosmology and astrophysics in LIM, and to compare the LIM
observables to galaxy detection for a number of surveys.
