Massive MIMO base stations use multiple spatial diversity branches, which are often assumed to be uncorrelated in theoretical work. Correlated branches are considered seldom since they are mathematically less tractable. For correlated Rician fading, only the first- and second-order moments have been explored. To describe propagation environments more accurately, full distribution functions are needed. This manuscript provides these functions for the maximum ratio combining effective channel, a quadratic form of a random complex normal channel vector. Its mean vector and covariance matrix are based on a plane wave model incorporating array geometry, antenna element pattern, power angular spectra and power delay profiles. Closed-form approximations of the distribution functions are presented, to allow the fast evaluation of many real-world scenarios. The statistical framework is used to show that low-directivity antenna elements provide better performance in angular constricted Rician fading with off-axis incidence than high-directivity elements. Moreover, two base station array layouts are compared, showing that a half-circle array illuminates a cell more evenly than a uniform linear array. With the full distribution functions available, performance can be compared over the full range of received powers and not only based on the average SNR.