Fringe fields in multipole magnets can have a variety of effects on the linear and nonlinear dynamics of particles moving along an accelerator beamline. An accurate model of an accelerator must include realistic models of the magnet fringe fields. Fringe fields for dipoles are well understood and can be modelled at an early stage of accelerator design in such codes as MAD8, MADX or ELEGANT. However, usually it is not until the final stages of a design project that it is possible to model fringe fields for quadrupoles or higher order multipoles. Even then, existing techniques rely on the use of a numerical field map, which will usually not be available until the magnet design is well developed. Substitutes for the full field map exist but these are typically based on expansions about the origin and rely heavily on the assumption that the beam travels more or less on axis throughout the beam line. In some types of machine (for example, a non-scaling FFAG such as EMMA) this is not a good assumption. In this paper, a method for calculating fringe fields based on analytical expressions is presented, which allows fringe field effects to be included at the start of an accelerator design project. The magnetostatic Maxwell equations are solved analytically and a solution that fits all orders of multipoles derived. Quadrupole fringe fields are considered in detail as these are the ones that give the strongest effects. Two examples of quadrupole fringe fields are presented. The first example is a magnet in the LHC inner triplet, which consists of a set of four quadrupoles providing the final focus to the beam, just before the interaction point. Quadrupoles in EMMA provide the second example. In both examples, the analytical expressions derived in this paper for quadrupole fringe fields provide a good approximation to the field maps obtained from a numerical magnet modelling code.