Arrays of Butler matrices provide a promising front-end design for massive MIMO transceivers with low cost and low complexity. However, this advanced design does not necessarily translate to effective applications, unless the angle-of-arrival (AoA) of signals avails to the Butler matrices. This paper presents an efficient approach to the unprecedented AoA estimation for the arrays of Butler matrices. Specifically, we design a new beam synthesis method to recursively narrow down and increasingly focus on the angular region of interest, and hence achieving robust estimation of the phase offset between Butler matrices. With the phase offset canceled in the received signals, we are able to identify the set of critical Butler beams with the dominating effect on the AoA estimation, and estimate the AoA accordingly with minimum signaling. The mean squared error of the proposed estimation is analyzed in the presence of non-negligible noises, with closed-form lower bounds derived. Validated by simulations, the proposed algorithm is able to indistinguishably approach the lower bounds, and significantly outperforms the state-of-the-art developed for discrete antenna arrays by orders of magnitude in terms of accuracy, especially in low signal-to-noise regimes.