We investigate the full counting statistics of a harmonically confined 1d short range Riesz gas consisting of N particles in equilibrium at finite temperature. The particles interact with each other through a repulsive power-law interaction with an exponent k > 1 which includes the Calogero–Moser model for k = 2. We examine the probability distribution of the number of particles in a finite domain [−W,W] called number distribution, denoted by N(W,N) . We analyze the probability distribution of N(W,N) and show that it exhibits a large deviation form for large N characterized by a speed N3k+2k+2 and by a large deviation function (LDF) of the fraction c=N(W,N)/N of the particles inside the domain and W. We show that the density profiles that create the large deviations display interesting shape transitions as one varies c and W. This is manifested by a third-order phase transition exhibited by the LDF that has discontinuous third derivatives. Monte–Carlo simulations based on Metropolis–Hashtings (MH) algorithm show good agreement with our analytical expressions for the corresponding density profiles. We find that the typical fluctuations of N(W,N) , obtained from our field theoretic calculations are Gaussian distributed with a variance that scales as Nνk , with νk=(2−k)/(2+k) . We also present some numerical findings on the mean and the variance. Furthermore, we adapt our formalism to study the index distribution (where the domain is semi-infinite (−∞,W]) , linear statistics (the variance), thermodynamic pressure and bulk modulus.