We propose a new method for decomposing electron density of a crystal into contributions associated with pair-wise chemical bonds. To this end, an ion-covalent crystal is represented using a neutral, closed shell cluster assembled from identical structural elements (SE) and embedded into the lattice electrostatic potential. The wave function of this cluster is calculated using the one determinant approximation. Then, a set of orthonormal, noncanonical, multicenter orbitals of the cluster valence states is generated, so as each orbital is localized on one structural element. The projection operators technique is used here, the valence molecular orbitals of the cluster being taken for the orthonormal basis set. In this construction, the first-order reduced density matrix of the cluster valence electrons is exactly the sum of the first-order reduced density matrices of the SE, and the latter is the exact sum of localized on this cluster orbitals densities. The localized orbitals are then transformed into directe