  We present a comprehensive method for determining stellar mass functions, and
apply it to samples in the local Universe. We combine the classical 1/Vmax
approach with STY, a parametric maximum likelihood method and SWML, a
non-parametric maximum likelihood technique. In the parametric approach, we are
assuming that the stellar mass function can be modelled by either a single or a
double Schechter function and we use a likelihood ratio test to determine which
model provides a better fit to the data. We discuss how the stellar mass
completeness as a function of z biases the three estimators and how it can
affect, especially the low mass end of the stellar mass function. We apply our
method to SDSS DR7 data in the redshift range from 0.02 to 0.06. We find that
the entire galaxy sample is best described by a double Schechter function with
the following parameters: $\log (M^{*}/M_\odot) = 10.79 \pm 0.01$, $\log
(\Phi^{*}_1/\mathrm{h^3\ Mpc^{-3}}) = -3.31 \pm 0.20$, $\alpha_1 = -1.69 \pm
0.10$, $\log (\Phi^{*}_2/\mathrm{h^3\ Mpc^{-3}}) = -2.01 \pm 0.28$ and
$\alpha_2 = -0.79 \pm 0.04$. We also use morphological classifications from
Galaxy Zoo and halo mass, overdensity, central/satellite, colour and sSFR
measurements to split the galaxy sample into over 130 subsamples. We determine
and present the stellar mass functions and the best fit Schechter function
parameters for each of these subsamples.
