We combine the 2MASS extended source catalogue and the 2dFGRS to produce an IR selected galaxy catalogue with 17,173 measured redshifts. We use this extensive dataset to estimate the J and K-band galaxy luminosity functions. The LFs are fairly well fit by Schechter functions with J: M*-5log h= -22.36+/-0.02, alpha= -0.93+/-0.04, Phi=0.0104+/-0.0005 h^3/Mpc^3 and K: M*-5log h= -23.44+/-0.03, alpha=-0.96+/-0.05, Phi=0.0108+/-0.0006 h^3/Mpc^3 (2MASS Kron magnitudes). These parameters assume a cosmological model with Omega=0.3 and Lambda=0.7. With datasets of this size, systematic rather than random errors are the dominant source of uncertainty in the determination of the LF. We carry out a careful investigation of possible systematic effects in our data. The surface brightness distribution of the sample shows no evidence that significant numbers of low surface brightness or compact galaxies are missed by the survey. We estimate the present-day distributions of B-K and J-K colours as a function of absolute magnitude and use models of the galaxy stellar populations, constrained by the observed optical and infrared colours, to infer the galaxy stellar mass function. Integrated over all galaxy masses, this yields a total mass fraction in stars (in units of the critical mass density) of Omega_*.h= (1.6+/-0.1)/10^3 for a Kennicutt IMF and Omega_*.h= (2.9+/-0.17)/10^3 for a Salpeter IMF. These values agree with those inferred from observational estimates of the star formation history of the universe provided that dust extinction corrections are modest.