The effects of finite particle number sampling on the net baryon number cumulants, extracted from fluid dynamical simulations, are studied. The commonly used finite particle number sampling procedure introduces an additional Poissonian (or multinomial if global baryon number conservation is enforced) contribution which increases the extracted moments of the baryon number distribution. If this procedure is applied to a fluctuating fluid dynamics framework one severely overestimates the actual cumulants. We show that the sampling of so called test-particles suppresses the additional contribution to the moments by at least one power of the number of test-particles. We demonstrate this method in a numerical fluid dynamics simulation that includes the effects of spinodal decomposition due to a first order phase transition. Furthermore, in the limit where anti-baryons can be ignored, we derive analytic formulas which capture exactly the effect of particle sampling on the baryon number cumulants. These formulas may be used to test the various numerical particle sampling algorithms.The goal of heavy ion collision experiments is to study the properties of very hot and dense matter. It has been proposed that in the most energetic collisions of nuclei at the RHIC and LHC a new state of matter, the so called Quark-Gluon-Plasma, has been created [1][2][3][4][5]. Lattice QCD results predict the transition at vanishing net baryon number density to be a crossover [6][7][8]. Due to the fermion sign problem, these calculations cannot be directly extended to the interesting region of high net baryon density where the cross-over may change to a first order phase transition. One of the main goals of current and future experimental programs is to find experimental signals for a possible first order phase transition and critical point in the phase diagram of the strong interaction, Quantum Chromo Dynamics (QCD).The systems created in these nuclear collisions are very small, rapidly expanding and not always in local thermal and chemical equilibrium. Thus the dynamical evolution of such a collision is far from trivial and sophisticated transport models are being employed. The current stateof-the-art of such models are the fluid dynamical hybrid models. In these models one uses a non-equilibrium initial state for a viscous fluid dynamical evolution, which is followed by a Boltzmann-transport description for the hadronic freeze-out phase. This setup is convenient as the fluid dynamical equations allow for a straight forward inclusion of the Equation of State (EoS). Since the systems created in heavy ion collisions are rather small (on the order of ∼ 10-20 fm) the application of standard fluid dynamics is limited. In fluid dynamics one usually assumes that the particle number in a given fluid element sufficiently large so that local fluctuations of the particle number, and hence of the baryon-and energy density, can be neglected. Since this is not strictly the case in heavy ion collisions new fluid * steinheimer@fias.uni-frankfurt.de † vko...