The formation of zonal flows from inhomogeneous drift-wave (DW) turbulence is often described using statistical theories derived within the quasilinear approximation. However, this approximation neglects wave-wave collisions. Hence, some important effects such as the Batchelor-Kraichnan inverse-energy cascade are not captured within this approach. Here we derive a wave kinetic equation that includes a DW collision operator in the presence of zonal flows. Our derivation makes use of the Weyl calculus, the quasinormal statistical closure, and the geometrical-optics approximation. The obtained model conserves both the total enstrophy and energy of the system. The derived DW collision operator breaks down at the Rayleigh-Kuo threshold. This threshold is missed by homogeneous-turbulence theory but expected from a full-wave quasilinear analysis. In the future, this theory might help better understand the interactions between drift waves and zonal flows, including the validity domain of the quasilinear approximation that is commonly used in literature. Singh et al. 2014). The WKE has the intuitive form of the Liouville equation for the DW action density J in the ray phase space (Parker 2016;Ruiz et al. 2016;Zhu et al. 2018c;Parker 2018; Zhu et al. 2018a,b):where Ω is the local DW frequency, Γ is a dissipation rate due to interactions with ZFs, and {·, ·} is the canonical Poisson bracket. (For the sake of clarity, terms related to external forcing and dissipation are omitted here.) However, as with all QL models, Eq. (1.1) neglects nonlinear wave-wave scattering, and in consequence, is not able to capture the Batchelor-Kraichnan inverse-energy cascade (Srinivasan & Young 2012) or produce the Kolmogorov-Zakharov spectra for DWs (Connaughton et al. 2015). Hence, a question remains as to whether the existing WKE for inhomogeneous turbulence can be complemented with a wave-wave collision operator C[J, J]. The goal of this work is to calculate C[J, J] explicitly. Starting from the generalized Hasegawa-Mima equations (gHME) (Krommes & Kim 2000; Smolyakov & Diamond 1999), we derive a WKE with a DW collision operator C[J, J] for inhomogeneous DW turbulence. Our derivation is based on the Weyl calculus (Weyl 1931), which makes our approach similar to that in Ruiz et al. (2016) where the QL approximation was used. However, in contrast to Ruiz et al. (2016), we do not rely on the QL approximation here but instead account for DW collisions perturbatively using the quasinormal approximation. The main result of this work are Eqs. (4.12) and (4.22). In this final result, DWs are modeled in a similar manner as in Eq. (1.1). The difference is that Eq. (4.12) includes nonlinear DW scattering, which is described by a wave-wave collision operator C[J, J] that is bilinear in the DW wave-action density J. The resulting model conserves the two nonlinear invariants of the gHME, which are the total enstrophy and the total energy. The present formulation is fundamentally different from the previously reported homogeneous weak-wave-turbulence mode...