The interlayer waveguide grating coupling efficiencies under angular (rotational) misalignments are simulated using the 3D rigorous coupled-wave analysis (3D-RCWA) together with the RCWA equivalent-index-slab (RCWA-EIS) method. As examples of conical diffraction, rotations about the two coordinate axes, x and z, defined by the vectors [1 0 0] and [0 0 1], respectively, as well as an arbitrary axis, defined by the vector [2 2 1], are simulated for binary rectangular-groove gratings. The interlayer grating coupling efficiency is approximated by the product of the top- and bottom-grating diffraction efficiencies (DEs). It is found that the bottom-grating DEs decrease about 25% when the bottom grating is rotated ±0.1 rad (5.73°) about the z-axis. DEs slightly increase (5% to 10% depending on the grating structures) when the bottom grating is rotated ±0.1 rad about the x-axis. This is consistent with the diffraction behavior of an over-modulated grating. When the bottom grating is rotated about the vector [2 2 1], the change in DEs is asymmetric with a 100% decrease at a rotation angle -0.1 rad and a 67% decrease at a rotation angle +0.1 rad. The method is shown to be computationally efficient and numerically stable for grating structures with optimized parameters, and the resulting bottom-grating diffraction efficiencies demonstrate similar trends as those calculated by the 3D finite-difference time-domain simulations. The procedure presented can be directly used in the analysis and design of interlayer waveguide grating coupling for optical interconnects in high-density integrated electronics.