AbstractQuantifying the contributions of distinct mineral populations in bulk magnetic experiments greatly enhances the analysis of environmental and rock magnetism studies. Here, we develop a new method of parametric unmixing of susceptibility components in hysteresis loops. Our approach is based on a modified Gamma‐Cauchy exponential model that accounts for variable skewness and kurtosis. The robustness of the model is tested with synthetic curves that examine the effects of noise, sampling, and proximity (similar coercivities) of susceptibility components. We provide a Python‐based script, the Hist‐unmix, which allows the user to adjust a forward model of up to three ferromagnetic components as well as a dia/paramagnetic contribution. Optimization of all the parameters is achieved through least squares fitting (Levenberg‐Marquardt method), with uncertainties of each inverted parameter calculated through a Monte Carlo error propagation approach. For each ferromagnetic component, it is possible to estimate the saturation magnetization (Ms), saturation remanent magnetization (Mrs) and the mean coercivity (). Finally, Hist‐unmix was applied to a set of weakly magnetic carbonate rocks from Brazil, which typically show distorted hysteresis loops (wasp‐waisted and potbellied). For these samples, we resolved two components with distinct coercivities. These results are corroborated by previous experimental data, showing that the lower branch of magnetic hysteresis can be modeled by the presented approach and might offer important mineralogical information for rock magnetic and paleomagnetic studies.