Photoionization of neutral molecules in the coma surrounding a comet produces heavy ions which contaminate and mass-load the solar wind. Cometary ion distribution functions in the vicinity of comet Giacobini-Zinner (G/Z) are calculated using a Monte Carlo method. The distribution function calculated behind the shock has both cold (≈2 keV) and hot (≈40 keV) components.