A new four-dimensional ocean interpolation system based on locally weighted least squares fitting is presented. A loess filter is used to interpolate irregularly spaced data onto a uniform grid. This involves projecting the data onto quadratic functions of latitude and longitude while simultaneously fitting annual and semiannual harmonics by weighted least squares. The smoothness scale of the mapping method adapts to match the data density, thus producing gridded estimates with maximum resolution. The filter has a vertical dimension, such that the data on adjacent levels are included in the computation. This greatly reduces the effects of discontinuities in data distributions between adjacent levels, since the estimates at each level are no longer independent. The loess scheme has been further modified so that the weighting of data points is adjusted to allow for the influence of both bathymetry and land barriers. This allows the bathymetry to influence the mapped fields in a natural way, reduces... Abstract A new four-dimensional ocean interpolation system based on locally weighted least squares fitting is presented. A loess filter is used to interpolate irregularly spaced data onto a uniform grid. This involves projecting the data onto quadratic functions of latitude and longitude while simultaneously fitting annual and semiannual harmonics by weighted least squares. The smoothness scale of the mapping method adapts to match the data density, thus producing gridded estimates with maximum resolution. The filter has a vertical dimension, such that the data on adjacent levels are included in the computation. This greatly reduces the effects of discontinuities in data distributions between adjacent levels, since the estimates at each level are no longer independent. The loess scheme has been further modified so that the weighting of data points is adjusted to allow for the influence of both bathymetry and land barriers. This allows the bathymetry to influence the mapped fields in a natural way, reduces...