A lattice theory is presented for homogeneous and heterogeneous systems containing molecules with orientation-dependent interactions. Correlations due to interactions are accounted for in first-order (quasi-chemical) approximation by allowing variations of the distribution of intermolecular contacts. To model heterogeneous systems, parallel lattice layers are allowed to be differently occupied. A partition function, written as a sum over distributions of molecules over orientations and locations and of intermolecular contacts, is derived for an arbitrary collection of monomeric species. Using a maximum term argument, self-consistent field equations are derived for the equilibrium distributions. These equations are solved numerically. Allowing lattice sites to be vacant, free-volume effects can be accounted for. Expressions are obtained for energy, entropy, chemical potentials, pressure, and surface tension. The capabilities of the method are illustrated by applying it to a number of specific systems, one of which exhibits a closed-loop coexistence curve. Another model in which one component is energetically anisotropic has an asymetric coexistence curve. Properties of coexisting phases and of the interface between them are investigated. With these models, it is found that orientational interactions may lead to a pronounced ordering in the interface and a negative interfacial entropy.