We present a computational method to calculate arbitrary pair correlation functions of an orthorhombic system in the most efficient way. The algorithm is demonstrated by the calculation of the radial distribution function of shock compressed liquid hydrogen.