A method based on Lucy's iterative algorithm is developed to invert the equation of stellar statistics for the Galactic bulge and is then applied to the K-band star counts from the Two-Micron Galactic Survey in a number of off-plane regions (10 deg>|b|>2 deg, |l|<15 deg). The top end of the K-band luminosity function is derived and the morphology of the stellar density function is fitted to triaxial ellipsoids, assuming a non-variable luminosity function within the bulge. The results, which have already been outlined by López-Corredoira et al., are shown in this paper with a full explanation of the steps of the inversion: the luminosity function shows a sharp decrease brighter than MK=-8.0mag when compared with the disc population; the bulge fits triaxial ellipsoids with the major axis in the Galactic plane at an angle with the line of sight to the Galactic centre of 12 deg in the first quadrant; the axial ratios are 1:0.54:0.33, and the distance of the Sun from the centre of the triaxial ellipsoid is 7860pc. The major-minor axial ratio of the ellipsoids is found not to be constant, the best fit to the gradient being Kz=(8.4+/-1.7)xexp(-t/(2000+/-920)pc), where t is the distance along the major axis of the ellipsoid in parsecs. However, the interpretation of this is controversial. An eccentricity of the true density-ellipsoid gradient and a population gradient are two possible explanations. The best fit for the stellar density, for 1300pc