The density field in the outskirts of dark matter halos is discontinuous due to a caustic formed by matter at its first apocenter after infall. In this paper, we present an algorithm to identify the splashback shell formed by these apocenters in individual simulated halos using only a single snapshot of the density field. We implement this algorithm in the code SHELLFISH (SHELL Finding In Spheroidal Halos) and demonstrate that the code identifies splashback shells correctly and measures their properties with an accuracy of $<5%$ for halos with more than 50,000 particles and mass accretion rates of $Gamma_textrm{DK14}>0.5$. Using SHELLFISH, we present the first estimates for several basic properties of individual splashback shells, such as radius, $R_textrm{sp}$, mass, and overdensity, and provide fits to the distribution of these quantities as functions of $Gamma_textrm{DK14}$, $ u_textrm{200m}$, and $z.$ We confirm previous findings that $R_textrm{sp}$ decreases with increasing $Gamma_textrm{DK14}$, but we show that independent of accretion rate, it also decreases with increasing $ u_textrm{200m}$. We also study the 3D structures of these shells and find that they generally have non-ellipsoidal oval shapes. We find that splashback radii estimated by SHELLFISH are $20%-30%$ larger than those estimated in previous studies from stacked density profiles at high accretion rates. We demonstrate that the latter are biased low due to the contribution of high-mass subhalos to these profiles and show that using the median instead of mean density in each radial bin mitigates the effect of substructure on density profiles and removes the bias.