cosmo surface dataset hdr

function ds = cosmo_surface_dataset(fn, varargin)