Convolution kernels for multiwavelength imaging
^{1} Institut d’Astrophysique Spatiale, CNRS, UMR 8617, Univ. ParisSud, Université ParisSaclay, IAS, Bât. 121, Univ. ParisSud, 91405 Orsay, France
email: alexandre.boucaud@ias.upsud.fr
^{2} Sorbonne Universités, UPMC Univ Paris 6 et CNRS, UMR 7095, Institut d’Astrophysique de Paris, 98bis bd Arago, 75014 Paris, France
^{3} Laboratoire des Signaux et Systèmes (Univ. ParisSud, CNRS, CentraleSupélec, Université ParisSaclay), 91192 GifsurYvette, France
Received: 9 June 2016
Accepted: 5 September 2016
Astrophysical images issued from different instruments and/or spectral bands often require to be processed together, either for fitting or comparison purposes. However each image is affected by an instrumental response, also known as pointspread function (PSF), that depends on the characteristics of the instrument as well as the wavelength and the observing strategy. Given the knowledge of the PSF in each band, a straightforward way of processing images is to homogenise them all to a target PSF using convolution kernels, so that they appear as if they had been acquired by the same instrument. We propose an algorithm that generates such PSFmatching kernels, based on Wiener filtering with a tunable regularisation parameter. This method ensures all anisotropic features in the PSFs to be taken into account. We compare our method to existing procedures using measured Herschel/PACS and SPIRE PSFs and simulated JWST/MIRI PSFs. Significant gains up to two orders of magnitude are obtained with respect to the use of kernels computed assuming Gaussian or circularised PSFs. A software to compute these kernels is available at https://github.com/aboucaud/pypher
Key words: methods: observational / techniques: image processing / telescopes / techniques: photometric
© ESO, 2016
1. Introduction
The pointspread function (PSF), also known as beam, is one of the main characteristics of any astronomical imager. It is a model of the diffraction pattern resulting from the interaction between the electromagnetic radiation and the instrument optics and detectors at every wavelength. Since most instruments operate on a single or a series of bandpasses (through e.g. filters), the resulting effective PSF is an integral of the monochromatic PSFs over the wavelength range, weighted by the instrumental throughput and the source energy distribution of a given astronomical object. A more accurate model can even include convolutional effects such as guiding errors, trailing effects from a scanning mode, smearing by the detector response, or even nonconvolutional effects like the brighterfatter effect. Once imaged, these model PSFs exhibit a complex shape, including anisotropy, wings, and spikes that extend far from the centre. Another classic feature of the PSF derived from the laws of optics is the radially oscillating pattern of the response, especially in the monochromatic case, creating a series of peaks and valleys. These secondary peaks can account for a nonnegligible amount of the total power of the PSF. For groundbased astronomy, however, the atmospheric turbulence creates a smearing that redistributes the power of these peaks and valleys and enables the PSF to be modelled by simple analytic profiles such as twodimensional (2D) Gaussians. On the contrary, space telescopes can benefit from a much higher resolution at the expense of a full complexity of the PSF. To cite a few examples, the effective PSF of IRAS maps was elliptical owing to the scanning strategy, so the angular resolution was strongly anisotropic, with ratios up to 1:6 (e.g. 0.75′ × 4.6′ at 25 μm, from Wheelock et al. 1994). More recently, the effective PSFs of the Planck/HFI^{1} maps appeared to have an ellipticity in the range 1.04 to 1.4, depending on the spectral band (Ade et al. 2011); the PSF of the PACS^{2} photometer (Poglitsch et al. 2010) on board the Herschel satellite, characterised by Lutz (2012), showed a narrow core, a trilobe pattern and knotty structured diffraction rings. As we push the boundaries of both optical performances and detector capabilities of future missions, optical designs highly increase in complexity. Hence, for upcoming space surveys (Euclid^{3}, WFIRST^{4}) or observatories (Athena^{5} or JWST^{6}), the characterisation and processing of elaborated PSFs become a crucial task.
Most astrophysical studies necessitate multiwavelength observations, either from multiple bands or filters within an instrument or from various instruments and telescopes. The different maps are affected by a particular PSF and the pixelbased data comparison cannot be straightforward. However, a technique widely used in multiband photometry is to perform the measurements on PSF homogenised data, that is to select a dataset as reference (usually the one with the worst resolution, or wider PSF) and transform the other datasets so they are PSFmatched with the reference PSF; this technique is called PSF homogenisation or PSF matching (see e.g. Bertin et al. 2002; Gordon et al. 2008; Darnell et al. 2009; Desai et al. 2012; Hildebrandt et al. 2012). Usually, PSF homogenisation is achieved by convolving the image with a kernel that is generated from the PSF corresponding to the image and the reference PSF. In the literature, one can distinguish between parametric kernels and nonparametric methods. Parametric methods use a fit of an analytic model to each PSF (Moffat, multiple Gaussians, etc.) or their decomposition on a proper basis (e.g. GaussHermite polynomials or shapelets), and result in an analytic expression of the kernel (e.g. Kuijken 2008; Hildebrandt et al. 2012; Liu et al. 2015); nonparametric methods use pixel information from the image (e.g. Alard 2000) or adopt effective PSF images (e.g. Gordon et al. 2008; Aniano et al. 2011) to compute the kernels.
With the purpose of taking the full complexity and angular extension of the PSFs of space instruments into account, we address the creation of PSFmatching kernels for multiwavelength studies. We then present two use cases for these kernels: one based on the Herschel satellite data and a second on simulations for the MidInfraRed Instrument Imager (MIRI) of the James Webb Space Telescope (JWST). We also deliver a programme called pypher that computes the kernels given two PSF images (see Appendix A). This code has initially been developed in preparation for the Euclid mission (Laureijs et al. 2010).
In Sect. 2, we describe the algorithm for the generation of convolution kernels used to match the resolution of images. In Sect. 3, we assess the improvement brought by these kernels on the multiwavelength study of dust properties with the Herschel satellite, and show in Sect. 4 the reconstruction power of such kernels on PSF simulations of the future JWST satellite, before summarising this work in Sect. 5.
2. Kernel generation
2.1. Data model
We first consider an astrophysical image y, observed with an instrument modelled as a linear invariant system, (1)where h is the PSF convolved with the unknown sky x, n is the image noise and ∗ stands for the discrete convolution (see e.g. Gonzalez & Woods 2008).
Given two PSF models h_{a} and h_{b}, where a and b refer to different frequency bands from the same or various instruments, the process we are interested in, referred to as PSFmatching, is to transform the image y_{a} acquired at the angular resolution of h_{a}(2)into an image y_{a,b} at the angular resolution of h_{b}where k_{a,b} is the matching kernel from h_{a} to h_{b}.
This paper presents a linear algorithm that computes the kernel k_{a,b} to produce the image y_{a,b} from the original image y_{a} through a convolution. To this end, we need to construct the kernel k_{a,b} such that (5)
2.2. Kernel generation
For such linear systems as Eq. (5), one can seek an estimate of k_{a,b}, denoted , which minimises the least squares criterion JHowever, the presence of the convolution makes the system illposed for the inversion, hence the solution to Eq. (7)is not stable. The only way to stabilise the solution is to add information. For the considered system, we use a technique called regularisation. We choose a ℓ_{2} norm to have a linear estimator and use Fourier filtering; and penalise the high frequencies in which we expect the noise to dominate, using a highpass filter d. This corresponds to adding a relative degree of smoothness between values of neighbouring pixels. (8)where d is the secondorder differential operator (i.e. 2D Laplacian matrix) (9)and μ the regularisation parameter, which tunes the balance between the data fidelity and the penalisation. Other norms, such as ℓ_{1}, ℓ_{2}ℓ_{1} or TV (Total Variations) are known to better preserve the image details, but produce non linear estimators that require iterative algorithms to solve.
Denoting the Fourier transform of any twodimensional vector u by , the convolution theorem states that the realspace convolution is equivalent to a termwise product in Fourier space (10)where ⊙ symbolises the termwise product between vectors/matrices.
Under these assumptions, the cancellation of the first gradient of the criterion (8)leads to the classical regularised mean square solution of Eq. (7)in Fourier space (11)where w is a Wiener filter with highfrequency penalisation (12)and stands for the complex conjugate of matrix .
The realspace convolution kernel k_{a,b}, is eventually obtained via the inverse Fourier transform of Eq. (11). For two instruments a and b, this kernel is thus only parametrised by the regularisation parameter μ. The optimal balance between the data and the penalisation is found by setting μ to the inverse of the signaltonoise ratio (S/N) of the homogenised image, y_{a} in this case.
We provide with this paper the pypher programme (see Appendix A); an implementation of the algorithm 1 presented below, which computes k_{a,b}. We note that the optical transfer function (OTF) that appears in the algorithm is the discrete Fourier transform of a signal that has been translated so that its peak value is the first vector entry (i.e. Im[0, 0] for an image).
3. Impact on Herschel data analysis
The pypher kernels allow us to convolve multiple astronomical images to a common angular resolution. In order to preserve all the information during this process, a good knowledge of the effective PSFs of the instruments used is required. We focus here on the particular analysis of Herschel photometric images that have been widely used in the last few years to measure the dust temperature and the spectral index β across many astronomical objects from multiband imaging with PACS and SPIRE.
Given the uncertainties on the PSFs of PACS and SPIRE instruments, since the beginning of the operational period of Herschel and in recent years, convolutions from PACS to SPIRE angular resolution have been performed assuming Gaussian PSFs with a given FWHM estimated from dedicated observations of asteroids. In 2011, a better characterisation of the PSF of the instruments allowed Aniano et al. (2011, ADGS11) to develop a method to construct convolution kernels assuming circular PSFs. While a full analysis of the PSF of the SPIRE instrument has been finalised recently by the Instrument Control Centre (ICC, Schultz 2015), a final PSF characterisation for the PACS instrument has not been released yet by the ICC (Lutz 2015). Parallel work has been performed by Bocchio et al. (2016), who computed PACS effective PSFs^{7} from the combination of Vesta and Mars dedicated observations, which will be used hereafter. The current knowledge of both PACS and SPIRE PSFs, and the pypher code allow us to construct effective convolution kernels with an unprecedented precision. In this section we show how and to what extent the use of different convolution kernels can affect the results.
3.1. Herschel PSFs and kernels
For each Herschel band, λ, we define four PSF images:

1.
The effective PSF (E_{λ}), taken from Bocchio et al. (2016) for PACS and Schultz (2015) for SPIRE.

2.
The Gaussian PSF (G_{λ}), computed fitting a 2D Gaussian profile to E_{λ}. The PSF G_{λ} has the same FWHM as E_{λ} but does not account for secondary lobes and faint structures.

3.
The circular PSF (C_{λ}), computed from E_{λ} by averaging the image intensity in annular bins. All information on the asymmetry of the PSF is then lost.

4.
The Aniano PSF (A_{λ}), a circular PSF from the ADGS11 paper.
For each of the first three types of PSF, we define the corresponding matching kernels between the bands λ_{1}and λ_{2} as (13)and compute them with pypher. We also consider the Aniano kernels defined here as and computed in ADGS11 via a different method than that used in this work.
The homogenised images can therefore be denoted by (14)
Fig. 1 Effective PACS 70 μm and SPIRE 350 μm PSFs on the left. The central column shows Gaussian, circular, and effective kernels for these PSFs computed following the procedure described in Sect. 3.1. The corresponding relative residual images (see Eq. (15)) are shown in the right column. The bottom kernel and residual image correspond to this work. All images are 120′′ × 120′′. 

Open with DEXTER 
3.2. Kernel comparison
In this paragraph we assess the impact of using approximations of effective PSFs in the homogenisation process, which is directly related to the creation of the convolution kernel. Because the effective kernel is the target of our kernel generation algorithm, it is expected to produce better results than the other types described above.
To compare the different types of kernels, we chose to measure the difference between the effective PSF at λ_{1} matched to λ_{2} and the effective PSF at λ_{2}. We define the relative residuals R^{X} as (15)for each type of kernel.
For this comparison test, we consider the matching of the PACS 70 μm PSF to the resolution of SPIRE 350 μm. In the remainder of this paragraph, we refer to these bands as 70 and 350.
The two images on the left of Fig. 1 represent the effective PSFs of PACS 70 μm E_{70} and SPIRE 350 μm E_{350}. The central column shows the kernels computed with pypher from Gaussian (top) and circular (middle) approximations of the effective PSFs, as described in Sect. 3.1, and those computed directly from the leftside PSFs (bottom). On these kernel images, one can see the characteristics of the input PSFs: the Gaussian approximation has a single lobe, the circular approximation is axisymmetrical and presents a second lobe, and the last approximation has two lobes and the general shape of E_{350}. Next to these kernels, the associated homogenisation residual images are shown (see Eq. (15)for computation).
Fig. 2 First (m(r), top panel) and second (σ(r), bottom panel) polar moments of the residuals (right column images of Fig. 1 + ) as a function of the distance r to the image centre. The different lines represent the kernel types used for the homogenisation: Gaussian (red dotted), circular (green long dashed), ADGS11 (black dashed) and effective (blue solid, this work). The vertical black dotted line indicates the HWHM of the SPIRE 350 μm PSF. 

Open with DEXTER 
Both and have residuals of the order of 10% within the first lobe of E_{350} (central region of radius equal to the half width at half maximum (HWHM) ≃ 12′′). Outside of that region, the reconstruction from these two kernels is even worse. In particular, the extinction ring that marks the transition between the first and second lobe does not exist in the Gaussian case and is slightly shifted owing to azimuthal averaging in the circular case, which leads to a big residual error in both cases (white circle on and ). Both kernel and residual images from ADGS11, which are not shown in Fig. 1, present very similar behaviour to the circular approximation. Using the kernel constructed with effective PSFs, , we obtain very homogeneous residuals of the order of 0.1%.
To analyse these residuals in more detail, we introduce the first m and second σ polar moments of the residuals as follows: (16)where R^{X}(r,θ_{i}) is the intensity of the residual image (computed using the kernel X) at a distance r from the image centre and at an angle θ_{i} = 2πi/n_{ang}, with n_{ang} = 100. The moments express the intensity and dispersion of the residuals along the PSF radius.
These two values, computed on the residual images for the four kernel types, are shown in Fig. 2. As previously stated, the Gaussian case (dotted red lines) is only stable within a circle of radius equal to the HWHM of E_{350}. At further distance it shows very high first and second moments close to unity, and establishes very poor matching. The circular case (longdashed green) exhibits a constant first moment below 10% at all distances, which is better than the Gaussian case. This is mainly due to the computation of the first moment that is very similar to the circularising process and averages out the measurements. However, the asymmetry and local structures of E_{70} and E_{350} are lost and the second moment is comparable to that of . To produce the kernels, ADGS11 used narrower versions of PACS and SPIRE PSFs than those used in this work. Both first and second moments (dashed black) therefore present bumps at the position of the lobes. Finally, the effective case (blue) best matches E_{350}. In amplitude, the first moment is ≲0.1% and the second moment ≲1% at all distances, even if we note that the second moment presents the two bumps observed earlier. By comparison, these moments are two orders of magnitude lower than those of the other three cases.
This concludes in a major improvement in using pypher kernels with effective PSFs with respect to Gaussian, circular or ADGS11 kernels from an image processing standpoint. Next we test again these kernels on the determination of meaningful parameters from data or simulations.
3.3. Dust properties study
Herschel observations are often used to retrieve information about dust properties in our Galaxy and in local galaxies. In this section we show how the choice of PSFs to construct the convolution kernels can affect pixelbypixel measurements.
We consider an edgeon galaxy at D ~ 10 Mpc from us with an intrinsic vertical profile given by (17)where z_{d} is the scale height of the vertical dust distribution. We introduce the angular distance θ ≃ z/D and its characteristic value θ_{d} ≃ z_{d}/D. Three scale heights are examined, z_{d} = 0.1, 1 and 10 kpc (scenarios a., b. and c. from Table 1), corresponding to θ_{d} ≃ 2,20 and 200′′. Their intrinsic dust profile is shown in Fig. 3.
Fig. 3 Intrinsic dust profiles (see Eq. (17)) for scenarios a., b. and c. listed in Table 1, corresponding to characteristic scale heights of z_{d} = 0.1, 1 and 10 kpc, respectively. These profiles are then convolved with Herschel PSFs to simulate an edgeon galaxy observed in different bands. 

Open with DEXTER 
Scale height, z_{d}, and signaltonoise ratio (S/N) adopted for the simulation of dust profiles of an edgeon galaxy.
We convolve the intrinsic dust profile to the effective PSFs of PACS 70, 100 and 160 μm and SPIRE 250 and 350 μm, while keeping the pixel size at 1′′. We then make the basic assumption that the dust in the whole galaxy has a temperature of T_{d} = 20 K and spectral index of β = 1.6 (following Planck Collaboration XI 2014), and rescale the convolved models accordingly.
In order to simulate real data, we add Gaussian statistical noise to the models and consider three different values of S/N, 10^{4}, 10^{2}, 10 (scenarios b., b.1 and b.2 from Table 1), with respect to the dust emission at the peak position (z = 0). These images are then homogenised to the resolution of SPIRE 350 μm using the four kernel types described in Sect. 3.1 and resampled with a common pixel size of 10′′.
Finally, using a minimum χ^{2} method, we fit the multiwavelength data on a pixelbypixel basis to a modified blackbody (18)where τ_{ν0} is the optical depth at the reference frequency ν_{0} and B_{ν}(T) is the blackbody radiation for a grain at temperature T.
Fig. 4 Parameter Δβ as a function of ΔT for the four kernel types, using an edgeon galaxy with the dust profile b. (z_{d} = 1 kpc, S/N = 10^{4}). For a given kernel type, there are 21 data points, each of which correspond to a measurement on a single pixel along the galaxy profile. The spread thus represents the systematic error on the β−T measurement for this kernel type. Error bars on both axes at the centre indicate the statistical errors on β and T obtained from the χ^{2}fitting routine. 

Open with DEXTER 
We define Δβ and ΔT as the deviations from the reference dust spectral index β and temperature T_{d}, respectively. Figure 4 shows Δβ as a function of ΔT for scenario b. where z_{d} = 1 kpc and S/N = 10^{4}, for the different kernel types. Except for the case where effective PSFs are used, systematic discrepancies are present (up to ΔT ± 1 K and Δβ ± 0.3), which are comparable in amplitude to the statistical errors (horizontal and vertical bars in Fig. 4), and a spurious strong negative correlation appears between dust temperature and β parameter.
In order to show where this effect is most significant within the galaxy profile, we illustrate in Fig. 5 the quantities ΔT and Δβ at various angular distances from the galaxy centre and for different kernel types (G: Gaussian; C: circular; A: ADGS11; E: effective). Each column thus represents a vertical cut of the modelled galaxy. Dashed lines indicate the characteristic scale height of the intrinsic dust abundance profile. The top three panels show scenarios a., b. and c. where the scale height varies and the S/N is kept constant at 10^{4}. The two bottom panels show scenarios b.1 and b.2 where the scale height is fixed at 1 kpc and we vary the S/N (see Table 1 for a summary). Each panel of Fig. 5 shows the relationship between Δβ and ΔT just as in Fig. 4, and adds the spatial information to the data points.
Depending on the considered dust scale height, deviations from the reference dust temperature and β parameter reach ±2 K and ±0.6, respectively, with higher deviations for shorter scale heights and at higher galactic latitudes. Regardless of the convolution kernel adopted, the negative correlation is clearly measured for all the pixels in the vertical cut.
Scenarios a., b. and c. show that very low (<1%) deviations from the reference values are obtained using the effective kernels, while the use of other convolution kernels leads to larger errors for z ≳ z_{d}. However, as expected, in decreasing the level of S/N (scenarios b.1 and b.2), the noise starts to dominate over the signal and very large discrepancies are observed in temperature and β, regardless of the adopted convolution kernel.
Fig. 5 Heat maps of ΔT and Δβ measured on a simulated edgeon galaxy at various angular distances θ from the centre, as a function of the kernel type used. Each panel represents a scenario from Table 1. Dashed lines indicate the scale height of the dust profile z_{d}. This figure highlights the very small scatter in both ΔT and Δβ when using the kernels with effective PSFs (column E on each plot, until the S/N becomes to low as in panel b.2). 

Open with DEXTER 
4. JWST PSF simulations
To show the reliability of our method for complexshaped PSFs, we now test our algorithm on simulated PSFs from a telescope with an uncommon optical design. For this purpose, we select the MidInfraRed Instrument Imager (MIRI) of the James Webb Space Telescope (JWST, Bouchet et al. 2015). As the JWST main mirror is made of several hexagonal mirrors, the optical response of this instrument shows highly nonsymmetrical features that need to be accounted for in the homogenisation process.
We used WebbPSF (Perrin et al. 2012), the official JWST PSF simulation tool to generate a set of four broadband PSFs of 5′′ × 5′′ centred at λ_{1} = 5.6 μm, λ_{2} = 11.3 μm, λ_{3} = 18.0 μm and λ_{4} = 25.5 μm in order to cover the spectral range of MIRI. These broadband PSFs were generated from a set of 20 monochromatic PSFs assuming flat spectral energy distribution of the source, and oversampled at four times the pixel size of the detector, corresponding to a pixel scale of 0.11 arcsec. They are shown in the first column (and top of the second column) of Fig. 6.
Using pypher, we then computed three matching kernels, namely, K_{λ1,λ4},K_{λ2,λ4} and K_{λ3,λ4} (following the notation from Sect. 3.1) to homogenise the first three PSFs to the angular resolution of E_{λ4}. Using the same procedure as in Sect. 3.2, we compared these homogenised PSFs to the original PSF using the residual formalism (15)applied to the MIRI bandpasses. The resulting residual images R_{λ1,λ4},R_{λ2,λ4} and R_{λ3,λ4} are shown on the last column of Fig. 6.
Fig. 6 Proof of concept of PSF homogenisation for the JWST/MIRI instrument. The PSFs from the first column at 5.6 μm, 11.3 μm and 18.8 μm, respectively, are homogenised to the PSF on top of the second column at 25.5 μm, using pypher kernels. The resulting homogenised PSFs are shown in the second column. They are visually indistinguishable from the effective one on top. The relative residuals on the right confirm that the central part of the PSF is reconstructed at least to 0.1%. 

Open with DEXTER 
The central region of the residual images, within the two main lobes of the E_{λ4} PSF, has a low level of residuals (~10^{5}) in the three configurations Along the image borders, there are some nonnegligible residual patches. A quick visual comparison with the homogenised PSFs (central column) shows that these patches correspond to extremely faint regions of the PSF (<10^{6} w.r.t. the peak) and thus have a very low impact in the matching process.
5. Conclusions
In this paper, we propose a new method for the generation of static PSF homogenisation kernels which is applicable for instruments presenting complex PSFs such as recent or future spaceborn telescopes. The PSF on such optical systems is hardly ever static over the fieldofview, but we restricted the purpose of this paper to the production of homogenization kernels for the study of regions of interest on the image, where the PSF can be considered nonvariable. The treatment of the PSF varying over the whole fieldofview of modern instruments cannot be linearized as in this work and requires a very different approach. It will be the subject of a following paper. The application on Herschel/PACS and SPIRE and JWST/MIRI instruments demonstrates the performance of the proposed algorithm in terms of low residuals (better than 10^{2}–10^{3} and 10^{5}–10^{6} for observed and simulated PSFs, respectively).
To assess the improvement brought by our algorithm for multiwavelength studies, we address the estimation of dust temperature and spectral index β of astronomical objects using multiband images taken in the submillimeter spectral range by Herschel. This estimation is made via pixelbypixel measurements across these images which have different intrinsic angular resolutions. Homogenisation kernels are thus traditionally used to bring all the images to the same angular resolution. Most of the analyses performed so far use either Gaussian kernels, or the circularised kernels produced by Aniano et al. (2011). However, effective PSFs of space imagers are anisotropic, so these methods are not accurate enough therefore to introduce systematic anticorrelation on β and temperature measurements with an amplitude that can be larger than the statistical noise. We have checked that using pypher kernels, systematic errors are in any case negligible compared to statistical noise.
Finally, we provide the pypher software (Boucaud 2016) to compute homogenisation kernels to be used for current and future instruments.
These PSFs are publicly available at http://idocherschel.ias.upsud.fr/sitools/clientuser/Herschel/projectindex.html.
Acknowledgments
A.B. would like to thank Hacheme Ayasso for useful discussions. We acknowledge the CNES (Centre National d’Études Spatiales) for supporting this work as part of the Euclid SGS (Science Ground Segment) within the Euclid Consortium. We acknowledge the Euclid Consortium, the European Space Agency and agencies and institutes supporting the development of Euclid. Part of this work has received funding from the European Union’s Seventh Framework Programme (FP7/2007–2013) for the DustPedia project (grant agreement No. FP7SPACE606847). This research has made use of the Herschel Database operated by the Integrated Data and Operation Center (IDOC) at the Institut d’Astrophysique Spatiale (IAS) under contract with CNES and CNRS. This research made use of Astropy, a communitydeveloped core Python package for Astronomy (Astropy Collaboration et al. 2013).
References
 Ade, P. A., Aghanim, N., Arnaud, M., et al. 2011, A&A, 536, A7 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Alard, C. 2000, A&AS, 144, 363 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Aniano, G., Draine, B. T., Gordon, K. D., & Sandstrom, K. 2011, PASP, 123, 1218 (ADGS11) [NASA ADS] [CrossRef] [Google Scholar]
 Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., et al. 2013, A&A, 558, A33 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Bertin, E., Mellier, Y., Radovich, M., et al. 2002, in Astronomical Data Analysis Software and Systems XI, 281, 228 [Google Scholar]
 Bocchio, M., Bianchi, S., & Abergel, A. 2016, A&A, 591, A117 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Boucaud, A. 2016, pypher: Python PSF Homogenization kERnels, http://dx.doi.org/10.5281/zenodo.61392 [Google Scholar]
 Bouchet, P., GarcíaMarín, M., Lagage, P.O., et al. 2015, PASP, 127, 612 [NASA ADS] [CrossRef] [Google Scholar]
 Darnell, T., Bertin, E., Gower, M., et al. 2009, in Astronomical Data Analysis Software and Systems XVIII, 411, 18 [Google Scholar]
 Desai, S., Armstrong, R., Mohr, J. J., et al. 2012, ApJ, 757, 83 [NASA ADS] [CrossRef] [Google Scholar]
 Gonzalez, R. C., & Woods, R. E. 2008, Digital Image Processing, 3rd edn. (Upper Saddle River, NJ: Prentice Hall) [Google Scholar]
 Gordon, K. D., Engelbracht, C. W., Rieke, G. H., et al. 2008, ApJ, 682, 336 [NASA ADS] [CrossRef] [Google Scholar]
 Hildebrandt, H., Erben, T., Kuijken, K., et al. 2012, MNRAS, 421, 2355 [NASA ADS] [CrossRef] [Google Scholar]
 Kuijken, K. 2008, A&A, 482, 1053 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Laureijs, R. J., Duvet, L., Sanz, I. E., et al. 2010, in SPIE Astronomical Telescopes+ Instrumentation, International Society for Optics and Photonics, 77311 [Google Scholar]
 Liu, J., Henning, C., Desai, S., et al. 2015, MNRAS, 449, 3370 [NASA ADS] [CrossRef] [Google Scholar]
 Lutz, D. 2012, PACS photometer point spread function [Google Scholar]
 Lutz, D. 2015, HerschelPACS document PICCMETN033, v2.2, http://herschel.esac.esa.int/twiki/pub/Public/PacsCalibrationWeb/bolopsf_22.pdf [Google Scholar]
 Perrin, M. D., Soummer, R., Elliott, E. M., Lallo, M. D., & Sivaramakrishnan, A. 2012, in SPIE Astronomical Telescopes+ Instrumentation, International Society for Optics and Photonics, 84423 [Google Scholar]
 Planck Collaboration XI. 2014, A&A, 571, A11 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Poglitsch, A., Waelkens, C., Geis, N., et al. 2010, A&A, 518, L2 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Schultz, B. 2015,Analysis Details of SPIRE Photmeter Beam Profiles, http://herschel.esac.esa.int/twiki/bin/view/Public/SpirePhotometerBeamProfileAnalysis2 [Google Scholar]
 Wheelock, S., Gautier, T., Chillemi, J., et al. 1994, IRAS Sky Survey Atlas ISSA, Explanatory Supplement, JPL [Google Scholar]
Appendix A: The pypher code
A Python code called pypher (Boucaud 2016), which computes the static PSF homogenisation kernels described in this work has been made publicly available and can be retrieved at: https://github.com/aboucaud/pypher
Once installed, this programme can be used through a commandline interface taking as input the PSF images (source and target) as fits files, and specifying the output filename for the kernel,
$ pypher psf_a.fits psf_b.fits kernel_a_to_b.fits
The tunable parameters are, first, the regularisation parameter μ of the Wiener filter (see Eq. (12)) that penalises the high frequencies, and should be set according to the image that will be homogenised, and second, the position angle of both PSFs with respect to their image to accurately take into account the PSF shape in the homogenisation process.
The programme takes less than a second on a single CPU to compute a kernel from two 512 × 512 PSF images.
All Tables
Scale height, z_{d}, and signaltonoise ratio (S/N) adopted for the simulation of dust profiles of an edgeon galaxy.
All Figures
Fig. 1 Effective PACS 70 μm and SPIRE 350 μm PSFs on the left. The central column shows Gaussian, circular, and effective kernels for these PSFs computed following the procedure described in Sect. 3.1. The corresponding relative residual images (see Eq. (15)) are shown in the right column. The bottom kernel and residual image correspond to this work. All images are 120′′ × 120′′. 

Open with DEXTER  
In the text 
Fig. 2 First (m(r), top panel) and second (σ(r), bottom panel) polar moments of the residuals (right column images of Fig. 1 + ) as a function of the distance r to the image centre. The different lines represent the kernel types used for the homogenisation: Gaussian (red dotted), circular (green long dashed), ADGS11 (black dashed) and effective (blue solid, this work). The vertical black dotted line indicates the HWHM of the SPIRE 350 μm PSF. 

Open with DEXTER  
In the text 
Fig. 3 Intrinsic dust profiles (see Eq. (17)) for scenarios a., b. and c. listed in Table 1, corresponding to characteristic scale heights of z_{d} = 0.1, 1 and 10 kpc, respectively. These profiles are then convolved with Herschel PSFs to simulate an edgeon galaxy observed in different bands. 

Open with DEXTER  
In the text 
Fig. 4 Parameter Δβ as a function of ΔT for the four kernel types, using an edgeon galaxy with the dust profile b. (z_{d} = 1 kpc, S/N = 10^{4}). For a given kernel type, there are 21 data points, each of which correspond to a measurement on a single pixel along the galaxy profile. The spread thus represents the systematic error on the β−T measurement for this kernel type. Error bars on both axes at the centre indicate the statistical errors on β and T obtained from the χ^{2}fitting routine. 

Open with DEXTER  
In the text 
Fig. 5 Heat maps of ΔT and Δβ measured on a simulated edgeon galaxy at various angular distances θ from the centre, as a function of the kernel type used. Each panel represents a scenario from Table 1. Dashed lines indicate the scale height of the dust profile z_{d}. This figure highlights the very small scatter in both ΔT and Δβ when using the kernels with effective PSFs (column E on each plot, until the S/N becomes to low as in panel b.2). 

Open with DEXTER  
In the text 
Fig. 6 Proof of concept of PSF homogenisation for the JWST/MIRI instrument. The PSFs from the first column at 5.6 μm, 11.3 μm and 18.8 μm, respectively, are homogenised to the PSF on top of the second column at 25.5 μm, using pypher kernels. The resulting homogenised PSFs are shown in the second column. They are visually indistinguishable from the effective one on top. The relative residuals on the right confirm that the central part of the PSF is reconstructed at least to 0.1%. 

Open with DEXTER  
In the text 