A novel kernel-based method is proposed for fast, highly accurate and numerically stable computations of polar harmonic transforms (PHT) in polar coordinates. Euler formula is used to derive a novel trigonometric formula where the later one is used in the kernel generation. The simplified radial and angular kernels are used in efficient computation PHTs. The proposed method removes the numerical approximation errors involved in conventional methods and provides highly accurate PHTs coefficients which results in highly improved image reconstruction capabilities. Numerical experiments are performed where the results are compared with those of the recent existing methods. In addition to the tremendous reduction in computational times, the obtained results of the proposed method clearly show a significant improvement in rotational invariance.