Source-linked AI summary
The Data Processing Pipeline for the MUSE Instrument
Peter M. Weilbacher, Ralf Palsa, Ole Streicher, Roland Bacon, Tanya Urrutia, Lutz Wisotzki, Simon Conseil, Bernd Husemann, Aurélien Jarno, Andreas Kelz, Arlette Pécontal-Rousset, Johan Richard, Martin M. Roth, Fernando Selman, Joël Vernet
TL;DR
MUSE’s complex 24-CCD observations and large raw datasets require dedicated automated reduction software. This paper describes the pipeline that calibrates raw exposures, converts them through pixel tables, and constructs datacubes, reporting broad scientific use while identifying calibration and variance-propagation caveats.
Problem
MUSE produces complex data from 24 CCDs and about 800 MB raw datasets, requiring a dedicated processing environment for automated reduction.
Method
The paper describes calibration and science-processing algorithms that transform raw CCD data into pixel tables and then fully reduced datacubes.
Results
More than 200 publications had produced new science results using MUSE data processed by the MUSE pipeline or ESO-provided datacubes.
Takeaways & Limitations
The pipeline’s single final resampling step supports variance propagation from raw data through the reconstructed datacube.
Takeaways & Limitations
Geometric calibration produces pseudo-pixels rather than exact sky sampling and must be corrected by astrometric calibration before final scientific analysis.
Abstract
from arXiv · showhide
Processing of raw data from modern astronomical instruments is nowadays often carried out using dedicated software, so-called "pipelines" which are largely run in automated operation. In this paper we describe the data reduction pipeline of the Multi Unit Spectroscopic Explorer (MUSE) integral field spectrograph operated at ESO's Paranal observatory. This spectrograph is a complex machine: it records data of 1152 separate spatial elements on detectors in its 24 integral field units. Efficiently handling such data requires sophisticated software, a high degree of automation and parallelization. We describe the algorithms of all processing steps that operate on calibrations and science data in detail, and explain how the raw science data gets transformed into calibrated datacubes. We finally check the quality of selected procedures and output data products, and demonstrate that the pipeline provides datacubes ready for scientific analysis.
1. Introduction
MUSE is a complex integral-field spectrograph whose large, multi-detector data volume motivated a dedicated, automated processing pipeline. The paper presents the pipeline’s design, processing algorithms, calibration requirements, and software implementation.
- MUSE is a medium-resolution, large-field integral-field spectrograph operated at ESO’s Very Large Telescope since October 2014.
- The pipeline was developed against requirements for image reconstruction, flux calibration, offsets, sky subtraction, exposure combination, and wavelength calibration.The stated targets included better than 1/4-pixel image precision, ±20% flux accuracy, and better than 1/20th-resolution-element wavelength accuracy.
- The instrument records data on 24 CCDs and supports wide-field and narrow-field observing modes, including adaptive-optics operation.Wide-field mode covers approximately 1′ × 1′, while narrow-field mode uses 25 mas pixels over about 7′.5 × 7′.5.
- Because raw data are complex and about 800 MB uncompressed per dataset, a dedicated processing environment was considered necessary.
- The paper describes the science-processing steps at a high level and then details the calibration and science procedures implemented in pipeline version 2.8.3.Version 2.8.3 was publicly released in June 2020 and supports all instrument modes, including NFM.
2. Science processing overview
MUSE science processing converts raw CCD exposures into pre-reduced pixel tables and then into fully reduced datacubes. The workflow combines CCD-level calibration, coordinate assignment, on-sky corrections, and cube reconstruction, with several optional steps.
- The pipeline splits processing into basic CCD-level calibration and post-processing that applies on-sky corrections and constructs the final datacube.Pre-reduced pixel tables connect the two processing levels.
- Raw science exposures are processed separately for all 24 CCDs, with data-quality and variance images initialized alongside each CCD image.Saturation and optional bad-pixel information are recorded before overscan and bias correction.
- Basic processing can include gain conversion, optional dark subtraction, wavelength cropping, illumination correction, and twilight-sky correction.The nominal useful wavelength range is 4750...9350 Å, or 4600...9350 Å in extended mode.
- Mandatory trace, wavelength, and geometry calibrations assign two spatial pseudo-pixel coordinates and one wavelength coordinate to each CCD pixel.This coordinate assignment produces a pixel table for each on-sky exposure.
- Post-processing merges IFU pixel tables, applies atmospheric-refraction and optional flux or telluric corrections, and performs sky subtraction before cube construction.The merged data are corrected using relative IFU efficiencies and a smoothed lamp-flat spectrum.
3. Full processing details
MUSE calibration begins with daytime master calibration data and is supplemented by twilight and nighttime on-sky calibrations. The detailed procedures characterize detector bias, noise, dark current, and pixel defects before science processing.
- MUSE calibration combines daytime calibration-unit exposures with twilight and nighttime on-sky calibrations.
- Master bias construction subtracts structured CCD bias patterns and combines individual bias frames using a default 3σ-clipped average.Bias levels include quadrant-dependent values and horizontal and vertical gradients.
- Read-out noise is estimated from statistics of 9 × 9-pixel boxes in bias-frame difference images and converted to electrons using the CCD gain.The measurements are made separately for each CCD quadrant.
- Master-bias application includes overscan correction, bias subtraction, and propagation of flagged pixels into the calibrated image.
- Dark exposures characterize dark current and identify hot pixels and bright columns; MUSE CCDs typically have about 1 e−h−1 pixel−1 dark current.A smooth dark model can be used to subtract pixel-level dark current without adding additional noise.
3.3. Flat-fielding
Flat-fielding corrects pixel-to-pixel sensitivity variations and identifies dark pixels, while tracing locates the 48 illuminated slices on each IFU CCD. The pipeline refines slice edges and fits smooth trace solutions from flat-field data.
- 3.3. Flat-fielding: Flat-field exposures are bias-subtracted and combined to correct sensitivity variations and locate dark pixels.Dark subtraction is optional and typically omitted because flat-field exposure times are short.
- 3.3. Flat-fielding: Dark pixels are identified as outliers within each traced slice, using thresholds of normally 5× the absolute median deviation below the median and 5σ above it.The upper threshold flags bright borders associated with electrons lost in dark columns.
- 3.3. Flat-fielding: Nominal-mode flat fields also flag weakly illuminated or negative pixels at the blue end of the CCD image.
- 3.4. Slice tracing: The tracing algorithm detects the 48 slice edges from averaged flat-field cuts, rejecting implausible widths between 72.2 and 82.2 pixels.If the expected slice count is not found, the edge threshold is iteratively adjusted.
- 3.4. Slice tracing: Gaussian fitting refines edge positions, after which rejected trace points are replaced by an iteratively clipped polynomial fit, typically of fifth order.The illustrated refinement shifted edge positions by less than 0.1 pixels.
3.5. Wavelength calibration
MUSE derives a separate two-dimensional wavelength solution for every slice from lamp-specific arc exposures, then corrects night-time data for temperature-driven wavelength shifts using sky emission lines.
- Computing the wavelength solution: A two-dimensional polynomial models each slice’s arc-line curvature and wavelength dispersion, with default horizontal and vertical orders of 2 and 6.The solution is evaluated at each CCD position to assign wavelengths.
- Arc-exposure preparation: Five exposures are typically taken for each HgCd, Ne, and Xe lamp, then combined lamp-wise with a 3σ-clipped average to reduce blends and improve faint-line detection.Separate lamp handling provides broader wavelength coverage while avoiding saturation of the brightest lines.
- Line identification: For each slice, lines are detected in a signal-to-noise spectrum from the central CCD column after median-background subtraction.The spectrum divides DATA by the square root of STAT, allowing detection whether or not arc frames were flat-fielded.
- Computing the wavelength solution: The final wavelength solution is obtained by iteratively fitting a two-dimensional polynomial to measured arc-line centers and reference wavelengths with 3σ clipping.Centroid uncertainties and scatter around the preceding one-dimensional fit can weight the solution.
- Night-time correction: Night-time observations can exhibit wavelength zero-point shifts of up to 1/10th of the spectral resolution because ambient temperatures differ from arc exposures.The pipeline estimates the shift from high-S/N sky spectra formed by averaging pixel-table values into 0.1 Å bins.
- Calibration context: The geometric calibration defines slice positions and the pipeline later applies astrometric calibration before producing the final science cube.Geometric calibration uses a pinhole mask to determine slice geometry, while astrometric calibration corrects global rotation and distortion.
3.7. Twilight sky-flat fielding
Twilight sky exposures provide a near-uniform illumination reference for correcting residual large-scale gradients and vignetting in MUSE data.
- Illumination reference: Bright sky-background exposures taken before sunset are used because their constant illumination closely resembles—but does not perfectly match—the illumination of objects on the sky.In adaptive-optics modes, notch filters block NaD laser light before further processing.
- Slice and IFU correction: The pipeline can use a raw illumination flat to correct relative illumination between slices within each IFU and save IFU sums for later relative scaling.The correction targets slice-to-slice illumination differences.
- Gradient correction: A smoothed, normalized white-light image is fitted with a two-dimensional model to isolate and correct large-scale gradients while suppressing slice-edge artifacts.The illuminated region is median-smoothed with a 5×7-pixel filter before fitting.
- Vignetting correction: For vignetting, residuals near the affected field edge are smoothed and modeled separately, using a default 4 × 4 polynomial or alternative Gaussian or median filtering.The procedure applies when a vignetting mask is supplied or Narrow Field Mode data are processed.
3.8. Line spread function measurement
The pipeline measures the line spread function from arc data to support sky subtraction, Raman-line modeling, and scientific analysis of reduced cubes.
- Purpose and scope: The LSF includes both the MUSE slice width and detector bin width because it is measured directly from CCD data.Accurate LSF characterization is important for sky subtraction and modeling laser-induced Raman lines.
- Inputs and preparation: Long arc-exposure sequences from all three lamps improve LSF sampling and signal-to-noise across the full wavelength range.The inputs and initial CCD corrections match those used for wavelength calibration.
- LSF models: The main implementation represents each IFU-and-slice LSF as an interpolated two-dimensional image varying with wavelength.An alternative represents the line shape with wavelength-dependent Gauss–Hermite coefficients stored in tables.
- Image-based fitting: For the image-based approach, isolated bright arc-line regions are placed in a special pixel table and fitted separately for each slice without combining images.Avoiding combination prevents biases during the LSF determination.
- Gauss–Hermite fitting: For the Gauss–Hermite approach, combined arc data are fitted in two stages: first widths and line fluxes, then the remaining coefficients with fluxes fixed.The two-stage fit separates line-strength optimization from line-shape modeling.
3.9. Sky subtraction
The MUSE pipeline supports sky subtraction from either the science exposure or an offset sky field, using emission-line and continuum modeling on largely unresampled data. It groups sky lines, fits their fluxes with an LSF-aware differential residual, and uses oversampled global sky spectra for subtraction.
- Sky-subtraction modes: MUSE supports sky subtraction from either the science exposure itself or a separately observed offset sky field.The algorithm is the same in both cases; only the creation and timing of the sky spectrum differ.
- Emission-line model: 7,800 OH-dominated sky lines are grouped into 52 wavelength groups, with lines sharing an upper level grouped because their fluxes approximately vary together.The input line list spans 3129–11000 Å and includes about 5,100 lines within the MUSE wavelength range.
- Emission-line model: The fit optimizes line-group flux factors and a linear wavelength correction by comparing modeled and measured neighboring-pixel intensity differences.The residual uses modeled and measured intensity differences, bin width, and estimated variance; each iteration includes the LSF and relevant sky-line fluxes.
- Sky-spectrum construction: The global sky spectrum is averaged over all sky regions, oversampled by four to 0.3125 Å bins, and modeled with a weighted average LSF.After the line fit converges, the continuum is obtained by subtracting the fitted lines from the averaged sky spectrum; users can override the computed continuum.
- Offset-sky processing: Offset sky fields are required when large objects fill the MUSE field because the sky continuum changes over tens of minutes and with moon illumination, while emission lines vary faster.The offset field is processed into a sky spectrum from selected unresampled pixels before fitting the line fluxes.
- Offset-sky processing: A 240 s offset-sky exposure yields estimated statistical S/N ≳250 for the continuum and >1,000 for moderately bright emission lines.The high S/N results from averaging the sky spectrum over a large part of the MUSE field.
3.10. Internal calibrations of science exposures
Internal science-exposure calibrations address Raman-scattered laser light and slice-level flux patterns. Raman correction models fitted molecular emission features, while self-calibration uses sky background to derive and apply segment-specific corrections before sky subtraction.
- Raman correction: Raman correction models laser-guide-star features from O2 near 6484 Å and N2 near 6827 Å, including faint secondary peaks extending about 50 Å.The features were identified after commissioning the adaptive-optics module in 2017 and are unresolved at MUSE’s spectral resolution.
- Raman correction: The Raman model assumes a constant ratio between each feature’s main and secondary peaks and uses molecular-physics line ratios with fitted absolute fluxes.Affected wavelength ranges are extracted, object and bad pixels are removed, and remaining sky pixels are fit using the LSF.
- Slice-level self-calibration: Flat-fielding leaves spatial structure at about 1–1.5% accuracy, including patterns from slice stacks and individual slices in wavelength-integrated images.The slice-level autocalibration is intended to address these residual illumination differences.
- Slice-level self-calibration: Self-calibration uses sky background to compute reference fluxes before sky subtraction, assuming an intrinsically flat background and sky in at least about one-third of each slice’s spaxels.The science exposure must therefore contain sufficient blank sky and no large objects covering the field.
- Slice-level self-calibration: The pipeline derives robust reference levels with repeated MAD-clipped and rejection-based averaging, then applies segment-specific correction factors to pixel data and variances.Segments are defined by slice and wavelength range, and the correction factors are applied in quadrature to pixel variances.
- Scope and controls: Deep-field self-calibration was tuned for HDF-S and UDF data, but other datasets can receive deviant segment factors that produce spectral discontinuities.Hand-optimizing the sky mask is identified as a way to improve results, and users can supply an external correction table.
3.11. Standard star handling
Standard-star processing derives throughput and telluric absorption corrections from reconstructed stellar cubes. The pipeline measures stellar fluxes, compares them with reference spectra, and produces smoothed response and telluric correction curves for each exposure.
- Calibration products: Spectrophotometric standard stars provide the instrument throughput and normalized telluric-absorption spectrum used to flux-calibrate science data.Throughput is represented by a response curve.
- Stellar processing: The standard-star recipe merges pixel tables, applies flat-field and atmospheric-refraction corrections, reconstructs a cube, and detects objects by thresholding.Detection thresholds range from 50 to 5σ until at least one object is found.
- Stellar processing: A Moffat point-spread-function fit is used for WFM stellar fields, with smoothing to suppress wavelength-to-wavelength noise and wavy continuum artifacts.Direct fits can otherwise produce unphysical differences between adjacent wavelength planes.
- Response and telluric curves: The measured selected-star flux is compared with an airmass- and extinction-adjusted reference spectrum to compute sensitivity as a function of wavelength.The sensitivity uses recorded counts, effective airmass, exposure time, and reference flux.
- Response and telluric curves: Telluric-affected wavelength ranges are interpolated across, and the resulting response curve is extrapolated, smoothed, and saved with telluric factors for each exposure.Piecewise cubic polynomials followed by a sliding average are usually more effective than a median filter for noise and outlier reduction.
3.12. Astrometric calibration
Astrometric calibration uses globular-cluster fields with Hubble Space Telescope catalogs to determine MUSE’s overall distortion and pixel scale. Matched stellar detections constrain a six-parameter world-coordinate solution with iterative rejection.
- Calibration fields: Astrometric fields use Milky Way globular clusters with existing Hubble Space Telescope imaging and catalogs containing a few hundred reference stars across a MUSE field.Outer cluster regions were selected for WFM, while denser central regions with a bright star were used for another mode.
- Astrometric solution: The calibration module applies standard on-sky preprocessing, reconstructs a cube, and optionally flux-calibrates the stellar field.Flux calibration is usually unnecessary when the astrometric input curves are unavailable.
- Astrometric solution: Typically about 100 matched objects constrain a six-parameter gnomonic-projection solution containing position, two scales, rotation, and shear.The fit is generally iterated twice with σ-clipping to reduce foreground-star effects.
- Pipeline controls: Object detection, centroiding, matching radius, iteration count, and rejection threshold are controlled through dedicated muse_astrometry parameters.Centroids may use Gaussian or Moffat fits.
3.13. Exposure offset calculation
The pipeline automatically measures relative offsets between MUSE exposures from detected sources, refines them iteratively, and writes offsets used during exposure combination.
- Offset measurement: Relative offsets are measured from reconstructed field-of-view images, using source lists detected independently in each exposure.The detection algorithm targets point sources and iteratively adjusts thresholds to reach predefined source-count limits.
- Iterative refinement: Pairwise source matches provide right-ascension and declination offsets that are refined through decreasing search radii.The initial estimate uses the mode of a two-dimensional offset histogram; later iterations select closer pairs and use their median offset.
- Global solution: A least-squares fit combines the measured pairwise field offsets into final offsets for each search radius.Optional weighting uses histogram peaks initially and the variance of selected relative offsets in later iterations.
- Application: The smallest-search-radius offsets are written to the offset list and applied to FITS-header field centers during exposure combination.Non-overlapping exposures receive zero relative offsets.
- Quality assessment: An approximate preview of the combined field of view is generated to support visual assessment of the computed offsets.An exposure map can also be created optionally.
4. Common algorithms
The pipeline’s common algorithms address resampling, cosmic-ray rejection, variance propagation, atmospheric refraction, and line-spread-function modeling. These procedures improve calibration and preserve useful uncertainty information while exposing important limits caused by resampling and modeling assumptions.
- Resampling and variances: A single final resampling step reconstructs the output cube, preserving per-voxel variance propagation while transferring variance into unstored covariances.The resulting voxel variances are accurate at the per-voxel level, but later aperture integrations are affected by covariance.
- Cosmic-ray rejection: 3D cosmic-ray rejection uses each pixel’s 26 neighbors and defaults to 15σ for single exposures or 10σ for multiple exposures.Tests found these cutoffs effective, with detected hits affecting approximately 0.06% of pixels in 600 s exposures and 0.11% in 1,500 s exposures.
- Resampling and variances: Drizzle is the default resampling method because it conserves object flux within the field of view.The input pixel scale is controlled with pixfrac independently in the three dimensions.
- Resampling and variances: Resampling creates spatially varying propagated variances, with values changing by up to a factor of ∼6 because of detector-to-cube geometry and nonlinear dispersion.Lower propagated variance corresponds to heavier local smoothing, while the patterns are geometric rather than changes in data trustworthiness.
- Resampling and variances: The resampled voxel-value distribution can appear non-Gaussian, but its wings decompose into nearly Gaussian components when voxels are grouped by similar propagated variance.These groups correspond to different amounts of local smoothing and reduced variance.
- Resampling and variances: Propagated variances describe the actual per-voxel noise level, but interpreting them as white noise underestimates aperture-based uncertainties because adjacent voxels are correlated.Directly measuring voxel-to-voxel fluctuations in a datacube cannot self-calibrate the variance level because those fluctuations are lower than the real noise.
- Resampling and variances: The variance estimator is limited when object photons contribute significantly to shot noise, requiring an object model with defined simulated noise in such fields.The approach is intended for cases where the objects of interest are faint enough for their photon contribution to be neglected.
- Atmospheric refraction: Atmospheric-refraction correction reduces a typical 1′′.2 blue-end shift at airmass 1.5 to below 0′′.04 with all four implemented refractive-index methods.The Filippenko formulation remains the default because it gives satisfactory results with MUSE data.
5. Implementation
The MUSE pipeline is implemented as ESO-compatible C modules built on CPL, with recipes usable through esorex, Reflex, and related systems. Its processing is parallelized across the 24 IFUs and post-processing loops, and the paper documents the complete reduction cascade and file flow.
- Software framework: The pipeline is written in pure C within the ESO software framework, using the Common Pipeline Library and supporting CFITSIO, wcslib, and FFTW.CPL provides the internal data structures and library foundation.
- User interfaces and integration: Reduction modules are shared-library plugins usable with esorex or Reflex, while esorex supports scripting and Reflex provides a graphical interface.The XML-based interface also supports documentation and integration with the MUSE-WISE processing system.
- Recipe architecture: Each reduction module, or recipe, calls a shared library to perform its task, with user-facing interfaces and basic C code generated from XML descriptions.The XML descriptions provide a common basis for documentation and system integration.
- Parallelization: Basic processing runs in parallel for the 24 IFUs, using up to 24 cores through OpenMP or externally launched processes.Post-processing is internally parallelized with OpenMP loops and can use all available cores.
- Distribution: The pipeline code is released as Open Source under GPL v2 through the ESO pipeline website.MUSE consortium releases are code-identical, although earlier versions used somewhat different packaging.
- Pipeline organization: The complete implementation cascade, including input and output file tags and all relevant recipes, is presented in Fig. 11.The figure includes both basic processing products and post-processing steps, including calibration modules.
6. Testing of the data quality
The pipeline was tested on simulated and publicly available data across calibration and science-processing tasks. These tests show accurate bias subtraction, wavelength and flux calibration, astrometry, sky subtraction, and near-ideal S/N scaling.
- Testing strategy: Testing used Instrument Numerical Model simulations during development and reproducible publicly available datasets for current performance checks.Simulations supported geometrical calibration, cosmic-ray rejection, and multi-exposure combination tests.
- Bias subtraction: Residual bias features are ≲0.1 e−h−1, below 1.5 × 10−21 erg s−1 cm−2 Å−1 for typical atmospheric transparency.A combined dark showed a vertical dark-current gradient of about 0.15 e−h−1 over 4000 pixels and a typical quadrant jump of about 0.05 e−h−1.
- Wavelength calibration: The wavelength calibration reaches below 0.024 Å, corresponding to ∼1 km s−1 velocity accuracy, and about 0.01 Å in the central wavelength range.The reported accuracy is preserved when resampling to typical 1.25 Å cube bins, while individual cube-line measurements are affected by sampling.
- Sky subtraction: Sky-subtraction residuals remain within 5% of the original sky level and generally within 1% in the continuum, except for limited line-spread-function cases reaching about 2%.The tests covered three publicly available exposures in different MUSE modes and examined the full wavelength range plus selected emission-line regions.
- Flux and astrometric calibration: Flux calibration is accurate to within 3–5%, while astrometric accuracy is about one-quarter of a spatial element or better.For 15 stars, the Gaia–MUSE separation was 0′′.036 ± 0′′.018, or ≲1/5 MUSE pixels; NFM Gaia validation was limited by insufficient star density.
- Exposure combination: After combining 55 exposures, S/N was 1.05 times the theoretical expectation in region p2 and 1.17 times it in region p1.S/N followed the ideal curve closely, although systematics were slightly stronger in region p1.
7. Conclusions & Outlook
The pipeline has supported more than 200 MUSE science publications and was developed through iterative prototyping, simulation, testing, collaboration, and complete software release. The authors identify several possible improvements and note applicability to future integral-field instruments.
- Scientific use: More than 200 publications had produced new science results using MUSE data processed by the pipeline or ESO-provided pipeline datacubes.Some studies additionally improved sky subtraction or background uniformity with tools such as ZAP, which was integrated into the ESO MUSE Reflex workflow.
- Development approach: The pipeline’s development benefited from early prototyping, raw-data simulation, hardware-test processing, sustained developer involvement, communication channels, and complete software publication.The simulation enabled a complete dry-run, while multiple institutes used the pipeline during hardware testing.
- Outlook: Potential improvements include nonlinear-detector corrections, improved telluric and sky-line treatment, line-spread-function propagation, and limited covariance propagation.These changes are presented as possibilities for improving data accuracy or subsequent analysis rather than as implemented pipeline capabilities.
- Outlook: Some algorithms described in the paper may also apply to future instruments such as ESO’s HARMONI and the possible BlueMUSE spectrograph.
Appendix A: Instrument modes
Appendix A lists the five MUSE observing modes and their nominal properties.
- Instrument modes: Table A.1 summarizes all five MUSE instrument modes and their nominal properties.
Appendix B: Instrument layout
Appendix B describes how MUSE distributes the focal-plane field across its integral-field units and slices before dispersion onto the detector.
- Instrument layout: The MUSE layout distributes the field among 24 integral-field units, each containing 48 slices.
- Instrument layout: Light from the slices is dispersed and recorded on the CCD, producing the step-pattern visible in the instrument data.The nominal element sizes differ between wide-field and narrow-field modes and vary slightly across the field.