The Monte Carlo method has been used to study a model system of 256 hard diatomic molecules, each consisting of two fused hard spheres of diameter
σ
with centres separated by reduced distance
L = L/σ
of 0.2, 0.4 and 0.6, at densities typical of the liquid state. The orientational structure of dense, hard diatomic fluids has been studied by calculating up to sixteen terms in the expansion of the total pair correlation function,
g
(
r
12
,
ω
1
,
ω
2
), in spherical harmonics. The coefficients
g
u'm
(
r
12) the series have been calculated as ensemble averages in the simulation. At short distances, the system exhibits a high degree of angular correlation, which increases with increasing density and elongation; however, this correlation is relatively short ranged at all densities and elongations, and in no case is there significant angular structure at distances greater than twice the major diameter of the molecule. In the nearest neighbour shell there is a strong preference for 'T-shaped’ pair orientations. At low elongations and densities the spherical harmonic coefficients are in close agreement with those predicted both by the ‘blip function’ theory and the solution of the Percus-Yevick equation for hard diatomics. The harmonic series for the total pair correlation function, is rapidly convergent at distances greater than
L
+
σ
, but slowly convergent at smaller distances. The results are suitable for use as a non-spherical reference system for perturbation calculations.