Point spread function (PSF) engineering has been pivotal in the remarkable progress made in high-resolution imaging in the last decades. However, the diversity in PSF structures attainable through existing engineering methods is limited. Here, we report universal PSF engineering, demonstrating a method to synthesize an arbitrary set of spatially varying 3D PSFs between the input and output volumes of a spatially incoherent diffractive processor composed of cascaded transmissive surfaces. We rigorously analyze the PSF engineering capabilities of such diffractive processors within the diffraction limit of light and provide numerical demonstrations of unique imaging capabilities, such as snapshot 3D multispectral imaging without involving any spectral filters, axial scanning or digital reconstruction steps, which is enabled by the spatial and spectral engineering of 3D PSFs. Our framework and analysis would be important for future advancements in computational imaging, sensing and diffractive processing of 3D optical information.