# spidast: a new modular software to process spectrointerferometric measurements

P. Cruzalèbes, Y. Rabbia, A. Jorissen, Alain Spang, S. Sacuto, E. Pasquato, A. Chiavassa, O. Chesneau, P. Fréville

► To cite this version:

P. Cruzalèbes, Y. Rabbia, A. Jorissen, Alain Spang, S. Sacuto, et al., spidast: a new modular software to process spectrointerferometric measurements. Monthly Notices of the Royal Astronomical Society, 2013, 432 (2), pp.1658-1671. 10.1093/mnras/stt593. hal-03526363

HAL Id: hal-03526363

https://hal.science/hal-03526363

Submitted on 14 Jan 2022

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# SPIDAST: a new modular software to process spectrointerferometric measurements*

P. Cruzalèbes,¹† Y. Rabbia,¹ A. Jorissen,² A. Spang,¹ S. Sacuto,³,⁴ E. Pasquato,² A. Chiavassa,¹,² O. Chesneau¹ and P. Fréville⁵

¹Laboratoire Lagrange, UMR 7293, Université de Nice-Sophia Antipolis, CNRS, Observatoire de la Côte d'Azur, Bd de l'Observatoire, B.P. 4229, F-06304 Nice cedex 4, France

²Institut d'Astronomie et d'Astrophysique, Université Libre de Bruxelles, Campus Plaine C.P. 226, Bd du Triomphe, B-1050 Bruxelles, Belgium

³Department of Astronomy and Space Physics, Uppsala Astronomical Observatory, Box 515, S-751 20 Uppsala, Sweden

⁴Institute of Astronomy, University of Vienna, Türkenschanzstraße 17, A-1180 Vienna, Austria

⁵Club d'Astronomie Copernic, Maison des Associations, 642 rue des Batteries, F-83600 Fréjus, France

Accepted 2013 April 4. Received 2013 April 4; in original form 2013 February 22

## ABSTRACT

Extracting stellar fundamental parameters from Spectrointerferometric (SPI) data requires reliable estimates of observables and with robust uncertainties (visibility, triple product, phase closure). A number of fine calibration procedures are necessary throughout the reduction process. Testing departures from centrosymmetry of brightness distributions is a useful complement. Developing a set of automatic routines called SPIDAST (made available to the community) to reduce, calibrate and interpret raw data sets of instantaneous spectrointerferograms at the spectral channel level, we complement (and in some respects improve) the ones contained in the AMDLIB Data Reduction Software. Our new software SPIDAST is designed to work in an automatic mode, free from subjective choices, while being versatile enough to suit various processing strategies. SPIDAST performs the following automated operations: weighting of non-aberrant SPI data (visibility, triple product), fine spectral calibration (subpixel level), accurate and robust determinations of stellar diameters for calibrator sources (and their uncertainties as well), correction for the degradations of the interferometer response in visibility and triple product, calculation of the centrosymmetry parameter from the calibrated triple product, fit of parametric chromatic models on SPI observables, to extract model parameters. SPIDAST is currently applied to the scientific study of 18 cool giant and supergiant stars, observed with the VLTI/AMBER facility at medium resolution in the K band. Because part of their calibrators have no diameter in the current catalogues, SPIDAST provides new determinations of the angular diameters of all calibrators. Comparison of SPIDAST final calibrated observables with AMDLIB determinations shows good agreement, under good and poor seeing conditions.

**Key words:** methods: data analysis – techniques: interferometric – stars: late-type.

## 1 INTRODUCTION

The power of optical–infrared interferometry to obtain information about the astronomical source morphology (including the angular size) is now well established. To properly determine the source properties, the quality of the measurements is an issue, which is still the subject of active research. As with any other measuring apparatus, the absolute calibration of the instrument (including

atmosphere) requires careful attention. Derived from the measurement of the mutual degree of coherence of the incident radiation field, on spatial frequencies sampled by the aperture-array configuration, the final interferometric observables are non-linear mixes of noisy quantities, and of parameter estimates with their own uncertainties.

In this paper, we propose to revisit and extend the existing data processing and calibration methods, in the aim to obtain reliable estimates and robust uncertainties for calibrated measurements of visibility and complex triple product. The careful reduction process, described in this paper, has been elaborated for the scientific study of a sample of 18 bright cool giant and supergiant stars (see Table 1). Measurements were obtained with the VLTI/AMBER facility at medium spectral resolution ($R = 1500$) in the K band, using triplets

*Applied on observations made with ESO telescopes at the Paranal Observatory under Belgian VISA Guaranteed Time programme ID 083.D-029(A/B), 084.D-0131(A/B) and 086.D-0067(A/B/C).

†E-mail: pierre.cruzalebes@oca.eu
<table><caption>**Table 1.** Science targets of our programme measured with VLTI/AMBER. Hipparcos parallaxes are in mas. The last column gives the calibrator(s) associated with each science target.</caption><thead><tr><th>Name</th><th>Spec. type</th><th>ω<sub>Hip</sub></th><th>m<sub>K</sub></th><th>Calibrator(s)</th></tr></thead><tbody><tr><td>α Car</td><td>F0II</td><td>10.6(6)</td><td>-1.3(3)</td><td>η Col/t Eri/HR 3282</td></tr><tr><td>β Cet</td><td>K0III</td><td>33.9(2)</td><td>-0.3(4)</td><td>η Cet</td></tr><tr><td>α TrA</td><td>K2II</td><td>8.4(2)</td><td>-1.2(1)</td><td>ε TrA/γ Lib/o Sgr</td></tr><tr><td>α Hya</td><td>K3II-III</td><td>18.1(2)</td><td>-1.1(2)</td><td>λ Hya</td></tr><tr><td>ζ Ara</td><td>K3III</td><td>6.7(2)</td><td>-0.6(2)</td><td>ε TrA/o Sgr</td></tr><tr><td>o<sub>1</sub> CMa</td><td>K2.5lab</td><td>0.2(4)</td><td>0.3(3)</td><td>HR 2411</td></tr><tr><td>δ Oph</td><td>M0.5III</td><td>19.1(2)</td><td>-1.2(2)</td><td>γ Lib/ε TrA</td></tr><tr><td>γ Hyi</td><td>M2III</td><td>15.2(1)</td><td>-1.0(4)</td><td>α Ret</td></tr><tr><td>o<sub>1</sub> Ori</td><td>M3III</td><td>5.0(23)</td><td>-0.7(2)</td><td>HR 2411</td></tr><tr><td>σ Lib</td><td>M3.5III</td><td>11.3(3)</td><td>-1.4(2)</td><td>51 Hya</td></tr><tr><td>γ Ret</td><td>M4III</td><td>7.0(1)</td><td>-0.5(3)</td><td>α Ret</td></tr><tr><td>L<sub>2</sub> Pup</td><td>M5IIIe</td><td>15.6(10)</td><td>-1.8(1)</td><td>HR 3282/η Col</td></tr><tr><td>CE Tau</td><td>M2lab-b</td><td>1.8(3)</td><td>-0.9(2)</td><td>40 Ori</td></tr><tr><td>T Cet</td><td>M5.5Ib/II</td><td>3.7(5)</td><td>-0.8(3)</td><td>γ Scl/t Eri</td></tr><tr><td>TX Psc</td><td>C7,2(N0)(Tc)</td><td>3.6(4)</td><td>-0.5(3)</td><td>θ Psc</td></tr><tr><td>R Scl</td><td>C6,5ea(Np)</td><td>2.1(15)</td><td>-0.1(1)</td><td>t Eri</td></tr><tr><td>W Ori</td><td>C5,4(N5)</td><td>2.6(10)</td><td>-0.5(4)</td><td>HR 2113/40 Ori</td></tr><tr><td>TW Oph</td><td>C5,5(Nb)</td><td>3.7(12)</td><td>0.5(4)</td><td>o Sgr/γ Lib</td></tr></tbody></table>

of 1.80-m auxiliary telescopes. Observations have been conducted during 15 observing nights between 2009 May and 2010 December (under Belgian VISA Guaranteed Time). The aim of this work is to test the data reduction and calibration software, on data with various qualities, that we started to develop in 2006.

## 2 DEFINING THE BASIC INTERFEROMETRIC OBSERVABLES

The coherent flux $c_{ij}$ is provided by the measurement of the instantaneous observables delivered by the AMBER data reduction software AMDLIB.<fnref n="1" /> This quantity traces the sine-like modulated component superimposed on the continuum component, in the observed intensity distribution. It is computed using a $\chi^2$ linear fit of the individual interferograms in the detector plane, for each spectral channel (see Tatulli et al. 2007 for details).

At this point, we define the other interferometric observables:

(i) the squared flux, which is the product of the photometric fluxes $p_i$ and $p_j$, associated with the baseline $B_{ij}$

$$f_{ij}^2 = p_i p_j; \qquad (1)$$

(ii) the squared visibility, which is the ratio of the squared modulus of the coherent flux to the squared flux

$$v_{ij}^2 = \frac{|c_{ij}|^2}{f_{ij}^2}; \qquad (2)$$

(iii) the complex bispectrum (Weigelt 1977), which is the product of the three coherent fluxes contributing to the baseline triplet ($B_{12}$, $B_{23}$, $B_{31}$)

$$b_{123} = c_{12}c_{23}c_{31}; \qquad (3)$$

(iv) the triple product, which is the ratio of the bispectrum to the product of the three fluxes contributing to the baseline triplet

$$t_{123} = \frac{b_{123}}{f_{12} f_{23} f_{31}}; \quad \text{and} \qquad (4)$$

(v) the closure phase, which is the argument of the complex bispectrum (or triple product)

$$\psi_{123} = \arg b_{123} = \arg t_{123}. \qquad (5)$$

## 3 REDUCING THE RAW DATA

To process the temporal series of frames (composing the 'exposures') produced by the AMBER instrument, we use the data reduction software AMDLIB, which consists of a core library of C functions plus a high-level interface in the form of a Yorick<fnref n="2" /> plugin (a high-level language environment, similar to IDL or MATLAB, which is in the public domain).

The standard procedure to extract the basic interferometric observables from the AMBER raw data includes a frame-selection step, which provides a set of 'good' frames according to a given criterion, or a sequence of several selection criteria (Malbet et al. 2011). The criterion used the most is based on the signal-to-noise ratio (SNR) of the fringe contrast in each frame. Two selection methods are provided in standard by AMDLIB:

(i) retaining a given percentage of frames sorted according to their SNR values (called the 'percentage' method) and

(ii) retaining the frames with SNR values higher than a given threshold (called the 'threshold' method).

With both methods, the choice of the optimal value (assumed not biasing the final results) of the selection criterion (in percentage of threshold) is done a posteriori from a sequence of processings using different values of the selection criterion (see Millour et al. 2007 for details). Two drawbacks of this are: the need for a several steps procedure (to be performed manually) and the risk for biasing the final results if the selection criterion is ill determined.

In our method, that we have implemented in a specific Yorick script added to AMDLIB, we apply the automatic procedure:

(i) remove only the aberrant measurements (in squared visibility and triple product), for each spectral channel independently. Using this specific script results in amounts of useful data larger than the amount obtained with the standard procedure, without increasing biasing effect (only the lower and upper tails of the input data histograms are rejected);

(ii) assign a weight to each individual measurement. After many preliminary tests, we decided to choose weights based on the SNR of the coherent flux and on the piston deviation; and

(iii) compute the temporal weighted average of the squared visibility and triple product for each spectral channel, with uncertainties derived from the bootstrap technique invented by Efron (1979).

These automatic operations mentioned above, needing neither several step selection nor manual choice of the optimal value of the selection criterion, are described in more details in the following subsections.

<footnote n="1">¹ www.mariotti.fr/data_processing_amber.htm</footnote>

<footnote n="2">² sourceforge.net/projects/yorick/</footnote>
## 3.1 Rejecting the aberrant data

To remove the aberrant measurements, we use a combination of physical and statistical basic criteria:

(i) strictly positive sum and product of the fluxes in the photometric channels, since negative flux measurements result from shortcomings in the determination of the continuum in each spectrum;

(ii) coherent-flux SNR higher than unity, since inclusion of data with low-SNR values reduces the reliability of the estimators of the interferometric observables (Millour et al. 2008); and

(iii) squared visibility greater than $Q1 - 1.5 \times \text{IQR}$ and lower than $Q3 + 1.5 \times \text{IQR}$ ('box-plot filtering'), where $Q1$ denotes the first quartile, $Q3$ the third quartile and $\text{IQR} = (Q3 - Q1)$ the interquartile range (Hoaglin, Mosteller & Tukey 1983). With triple product data, this third criterion applies to the squared visibility for each baseline, and box-plot filtering on $\tan \psi_{123}$ is added to the list of criteria.

Based on the removal of aberrant measurements at the spectral channel level, rather than on the selection of the 'best' frames, our approach provides an amount of useful data larger than the amount obtained with the standard frame-selection procedure. For data shown in Fig. 1, obtained for a 'good' seeing on the calibrator target $\gamma$ Lib, the number of removed input data varies only slightly from one spectral channel to the other (typically less than 3 per cent). Besides, various amount of data are rejected, according to the seeing conditions.

Using the seeing parameters $\epsilon_0$ (seeing angle) and $\tau_0$ (coherence time), the percentages of rejected data, respective to three baselines

<figure><img src="image_3.png" /><figcaption>**Figure 1.** Spectral distribution (in K) of the number of non-aberrant data, remaining after box-plot filtering, from an exposure of 300 input frames with the calibrator $\gamma$ Lib. From top to bottom: baselines A0–D0 (32 m), D0–H0 (64 m), A0–H0 (96 m) and triplet A0–D0–H0. Blue lines: good seeing conditions (seeing angle = 0.5 arcsec, coherence time = 13 ms). Red lines: poor seeing (1.7 arcsec, 1.3 ms).</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>**Figure 2.** Squared visibility at 2.2 µm obtained with AMDLIB (black solid lines + error bars) w.r.t. the percentage of selected frames sorted according to their SNR value, for the calibrator $\gamma$ Lib (baseline A0–D0). Top panel: good seeing conditions. Bottom panel: bad seeing. Red lines: values obtained with our method (solid line: central value; dashed lines: upper and lower limits given by the associated uncertainty).</figcaption></figure>

and two seeing conditions are given in the following lines (where 'good seeing' refers to $\epsilon_0 = 0.5$ arcsec, $\tau_0 = 13$ ms, and 'poor seeing' refers to $\epsilon_0 = 1.7$ arcsec, $\tau_0 = 1.3$ ms):

(i) shortest baseline A0–D0 (32 m): 2 per cent (good seeing) to 9 per cent (poor seeing);

(ii) intermediate baseline D0–H0 (64 m): 2 per cent to 16 per cent;

(iii) longest baseline A0–H0 (96 m): 2 per cent to 25 per cent; and

(iv) in addition, for the baseline triplet A0–D0–H0, the rejection of data amounts 7 and 50 per cent.

Fig. 2 shows the variation of the final squared visibility (at 2.2 µm) produced with AMDLIB, w.r.t. the percentage of selected frames, based on the SNR and pertaining to the calibrator target $\gamma$ Lib, observed with the shortest baseline A0–D0 (32 m) in good- and poor seeing conditions (top panel and bottom panel, respectively). The squared visibility produced with our method is shown for comparison.

In this example, the position of the plateau of each curve indicates a value of the selection threshold of 50 per cent, giving the a posteriori determination with AMDLIB. There we find that our automatic method gives a squared-visibility value close to the 50 per cent value, without a need for the several steps procedure, in good and poor seeing conditions. Besides, Fig. 3 shows that the standard AMDLIB frame-selection method keeps the extremal $V^2$ values in poor seeing condition, which may bias the result of the calculation of the average over the frames contained in each exposure. On the contrary, the box-plot filtering used with our method rejects only the
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Histograms of the squared visibilities at 2.2 µm, obtained with the calibrator $\gamma$ Lib (baseline A0–D0), under good (top panel) and poor (bottom panel) seeing conditions. Filled turquoise bars: final histogram after aberrant-data rejection (our method). Full black steps: initial histogram (no selection). Dashed red steps: final histogram after 50 per cent frame selection with AMDLIB.</figcaption></figure>

lower and upper tails of the $V^2$ histogram, thus, leading to estimates of the temporal average, more reliable than with the frame-selection method.

Next subsection describes the weighting scheme used to average the non-aberrant data over each exposure. Note that selecting the frames with the standard AMDLIB procedure is equivalent to assigning unity weights to the selected frames and null weights to the others.

<figure><img src="image_6.png" /><figcaption>**Figure 4.** Interferometric observables given by our reduction code with the calibrator $\gamma$ Lib, before calibration (good seeing). Top-left panel: squared visibilities for the baselines A0–D0 (blue), D0–H0 (green) and A0–H0 (red) w.r.t. the channel index in K. Top-right panel: real (blue) and imaginary (red) parts of the triple product for the baseline triplet A0–D0–H0. Because of their low levels, uncertainties in the bottom panels are drawn using histograms, rather than error bars (with the same colour code as the top panels).</figcaption></figure>

## 3.2 Computing the raw observables

To reduce the instrumental effects on the data, *at the spectral channel level*, we compute in each exposure (300 frames in our example) the weighted average of the squared visibility, and of the complex triple product. For the visibility, as weight associated with each frame we use the ratio of the SNR to the relative excursion of the piston, where the relative excursion of the piston itself is the ratio of the piston excursion (absolute value) for a given frame, to the average over the whole exposure. For the triple product, each individual weight is given by the geometrical mean of the three weights associated with the three baselines. The final uncertainties are derived from the variance of the distribution of the weighted means, obtained by random sampling with replacement of the original series of data (direct bootstrap method) (Efron & Tibshirani 1993). Note that our specific script computes the real and imaginary parts of the triple product, from which we extract the calibrated closure phase (see Section 4.2).

Fig. 4 shows an example of squared visibility and triple product produced by the reduction script (before calibration) from one single exposure of 300 input frames, obtained with the calibrator target $\gamma$ Lib, for the three baselines A0–D0 (32 m), D0–H0 (64 m) and A0–H0 (96 m).

In complement to these interferometric observables, we compute two other quantities used further in the calibration procedure (see Section 4.2.3):

(i) the weight of each exposure defined as the product of two ratios: one is the average of the SNR to its variance; the other is the average of the inverse of the piston excursion (absolute value) to its variance and

(ii) the final number of non-aberrant data, for each spectral channel.

<figure><img src="image_7.png" /></figure>
# 4 CALIBRATING THE DATA

The calibration process is the key point to obtain accurate estimates of the ‘true’ (intrinsic) observables of the science targets. To correct the measurements for environmental and instrumental instabilities, we observe reference stars (the calibrators), with known angular diameters and independently well-defined brightness distributions.

To calibrate the wavelength-dependent measurements obtained with VLTI/AMBER, we have developed a library of iDL functions, included in the Spectrointerferometric Data Analysis Software Tool (SPIDAST) modular software suite (Cruzalèbes, Spang & Sacuto 2008; Cruzalèbes et al. 2010), which allows one to:

(i) link each spectral channel to a wavelength value. We use a specific method based on the correlation of calibrator measured spectra with synthetic templates given by the MARCS model (Gustafsson et al. 2008) and

(ii) measure and correct for the degradation of the spatial coherence. As calibrators, we use stars with brightness distribution of limb-darkened (LD) discs, for which angular diameters are given by fits of synthetic spectra on the wide-band spectrophotometric measurements in the infrared (Cruzalèbes et al. 2010).

## 4.1 Computing the spectral shifts

To reach a high precision in the angular diameter estimation using model fitting, spectral calibration is a critical point (see e.g., Domiciano de Souza et al. 2008; Wittkowski et al. 2008; Stefl et al. 2011). Given the lack of any internal instrumental module for wavelength calibration in the optical setup of AMBER (Petrov et al. 2007; Robbe-Dubois et al. 2007), AMDLIB provides calibrated wavelength tables, computed from a theoretical polynomial dispersion law (Mérand et al. 2010), but with coefficients still badly known, leading to wavelength shifts up to ~10 pixels in the detector plane, at medium spectral resolution (Malbet et al. 2011). To correct for this drawback, AMDLIB shifts the wavelength table using the correlation of the measured spectrum with a template table containing the telluric lines.

To improve the precision of the wavelength calibration, observers usually compute the coefficients of the dispersion law by identifying visually some prominent telluric lines in the measured spectrum (e.g. Kraus et al. 2009; Ohnaka et al. 2011; Weigelt et al. 2011). We propose an alternative automatic approach, based on the cross-correlation of the calibrator measured spectra with synthetic templates produced by the MARCS + TURBOSPECTRUM codes (Alvarez & Plez 1998; Plez 2012), including the telluric lines as well. This method was previously suggested by Plez (2003) for Gaia, and by Decin et al. (2004) or Decin & Eriksson (2007) to calibrate the Spitzer spectrograms. It is particularly suitable for stellar spectra showing easily identifiable spectral features, as it is the case with our sample of cool giant calibrators, showing strong CO bands in the 2.126–2.474 µm spectral range (Martí-Vidal et al. 2011). Fig. 5 shows the synthetic spectrum produced by the MARCS + TURBOSPECTRUM code, for the calibrator γ Lib. Some reference spectral lines are shown for identification of the corresponding stellar atmospheric elements.

For each exposure, obtained with each calibrator, and with the working instrumental setting, our automatic spectral calibration process contains the following steps:

(i) apply the heliocentric and systemic radial velocity corrections, and multiply the synthetic spectrum by the atmospheric and instrumental transmittance profiles. Atmospheric transmission data

<figure><img src="image_8.png" /><figcaption>**Figure 5.** Spectral distribution of the synthetic model surface flux of γ Lib, produced by the MARCS + TURBOSPECTRUM code in K. Dash-dotted line: blackbody continuous spectrum with $T_{eff} = 4660$ K. Dash-triple dotted line: Engelke continuous spectrum with the same effective temperature.</figcaption></figure>

for the southern sites (CTIO, Chile) are given by the USAF atmospheric code PLEXUS (Cohen, Wheaton & Megeath 2003);

(ii) remove the continuum parts of the raw and synthetic spectra, and normalize the resulting spectra. As a first approximation, we estimate the continuum by decreasing the spectral resolution to $R = 40$. The main drawback of this method is to produce apparent ‘pseudo-continua’ that are lower than the real continuum levels, even in almost line-free regions (Rix et al. 2004). Since our goal is only to perform spectral calibration by correlation, this drawback has no effect on the final result;

(iii) divide the spectrum in contiguous subwindows of the same size, and compute the wavelength shift by phase correlation (Vera & Torres 2008) in each subwindow. If the spectrum is divided in *n* subwindows, we apply a polynomial law of degree *n* − 1 to correct for the wavelength shifts. Subpixel precision is obtained by embedding the cross-power spectrum in the middle of a two times larger array of zeros, before computing the phase-correlation function, by inverse Fourier transform of the final cross-power spectrum (Guizar-Sicairos, Thurman & Fienup 2008);

(iv) finds the position of the barycentre of the correlation peak for each subwindow, computes the coefficients of the interpolation law (assumed to be polynomial) and corrects the global wavelength table; and

(v) calculates the residuals of the wavelength shifts, using once again the phase-correlation method on the wavelength-corrected spectrum.

For each observing night, the final wavelength table, associated with the working instrumental setting, is obtained by ensemble average of the corrected tables (of all exposures with all calibrators measured with the same instrumental setting), rejecting the tables showing wavelength-shift residuals higher than 8 nm (~5 times the nominal spectral resolution). Fig. 6 shows the final continuum-corrected spectrum, given by the spectral-calibration procedure, with the calibrator γ Lib. The model spectrum drawn with the dashed line is produced by the MARCS + TURBOSPECTRUM code.

To show the difference between the wavelength calibration provided by SPIDAST and AMDLIB, we plot in Fig. 7 the corrections from the initial polynomial dispersion law of AMBER, given by the two software packages. While AMDLIB computes only a unique correction value applied over the whole spectrum, our method computes
<figure><img src="image_9.png" /><figcaption>Figure 6. Spectral distribution of the continuum-subtracted flux of γ Lib measured with VLTI/AMBER, after spectral calibration. Dashed line: model spectrum (including tellurics).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 7. Wavelength correction from the initial polynomial dispersion law of AMBER, obtained with SPIDAST (in red) and with AMDLIB (in blue), for each spectral channel of the K band (first observing night).</figcaption></figure>

a polynomial correction (here parabolic from three subwindows), with subpixel deviations with respect to the polynomial law, caused by the final ensemble averaging of the corrected wavelength tables. The choice of the degree of the polynomial correction provided by SPIDAST allows the accuracy of the wavelength correction to be adapted to the scientific programme, reaching a subpixel precision if necessary.

## 4.2 Calibrating spectrointerferometric data

In order to assess and correct for the measurement defects on the science targets, proper interferometric calibration (not yet fully supported with AMDLIB) is needed. The Instrument Response Function (IRF), also called ‘transfer function’ (Perrin 2003; van Belle & van Belle 2005; Boden 2007; Cruzalèbes et al. 2010), is derived from observations of calibrators as concomitant as possible to the ones of the science targets.

The standard calibration method is based on the assumption that a given ‘raw’ (measured) quantity q of any observable Q, is equal to the product of the ‘true’ (intrinsic) quantity Q, multiplied by a global degradation factor R_Q (in other words the IRF). Such a factor cannot practically be modellized with sufficient reliability (numerous parameters and some of them remaining unknown).

To calibrate the quantity related to the science target, a second assumption applies, as soon as appropriate conditions (described

later on) are satisfied and can be then determined by using the relation

$$ \frac{q_{\text{sci}}}{Q_{\text{sci}}} = \frac{q_{\text{cal}}}{Q_{\text{cal}}} = R_Q, \qquad (6) $$

with self-explanatory notations, stating that the degradation factor is identical for both the calibration source and the science target. This assertion is all the more valid that the required conditions are fulfilled. The relation here above yields an empirical determination of Q_{sci}, via R_Q, relying on commonplace statistical tools. The main condition to satisfy is that the IRF is stable enough, and the standard approach to avoid too large variations of the IRF is to perform sequences, where observation of calibrator and science targets are interlaced and repeated several times, rapidly sampling the IRF. So, the actual stability requirement is that the IRF might be slowly variable over a given ‘calibrator-science-calibrator’ observation sequence. Besides, the calibrator and the science target should be close enough regarding time of observation, angular separation and brightness (in the spectral domain of work), so that adjustments of the whole interferometer would not be significantly modified from one source to the other. In these conditions, only the atmospheric turbulence might cause substantial degradations of the measurements. For this latter case, servo-loop systems are to be used to reduce turbulence effects.

Usually, calibrators are found in dedicated catalogues and the corresponding sources should meet the desired requirements. However, the number of calibrators is necessarily limited and somewhat depends on the type of science targets, so that the selected calibrators might stand out from the ideal conditions: point-like source, small angular separation and brightness matching.

In practice, the limited choice of calibrators makes it necessary to accept ‘partially’ resolved targets (van Belle & van Belle 2005), angular separations counted in degrees, and brightness mismatch amounting some units of magnitudes. For example, a bright calibrator is rarely found in a close angular vicinity of a given bright science target. Such a situation nevertheless remains convenient for the calibration procedure, since switching between science targets and calibrators is fast enough and does not affect the configuration of the interferometer, and its response as well.

This procedure starts with a first selection of calibrators via the SearchCal tool,<fnref n="3" /> created by the JMMC working group ‘Catalogue of calibration sources’, providing rough estimates for angular diameters, with comparatively large uncertainties. Following this first selection, the sources are searched in the calibrator catalogues of Bordé et al. (2002) and Mérand, Bordé & Coudé Du Foresto (2005), in order to find angular diameters of better quality. Since some calibrators do not have a diameter in the published catalogues, it is necessary to make our own determinations of those diameters (see Section 4.2.1). Moreover, in order to control these determinations and to build a homogeneous set of diameters, we apply our own method to recalculate all the diameters for the selected calibrators. The good agreement found between our determinations and the ones of the Bordé’s Catalogue of Calibrator Stars for Long-Baseline Stellar Interferometry (CCSL; as shown in Table 2 and in Fig. 8) attests to the reliability of our specific determinations of angular diameters.

### 4.2.1 Determining calibrator angular diameters

To derive the IRF for each spectrointerferometric (SPI) observable, we first need to get a reliable estimate of the angular diameter of

<footnote n="3">³ www.mariotti.fr/searchcal_page.htm</footnote>
<table><caption>**Table 2.** Angular diameter values (in mas) of the suitable calibrators, found in the CCSL, derived from the final fit of the MARCS model spectrum and found in the JSDC. The third column gives the effective temperature (in K) deduced from the spectral type, used for the model.</caption><thead><tr><th>Name</th><th>Spec. type</th><th><i>T</i><sub>eff</sub></th><th><i>Φ</i><sub>CCSL</sub></th><th><i>Φ</i><sub>final</sub></th><th><i>Φ</i><sub>JSDC</sub></th></tr></thead><tbody><tr><td>α Ret</td><td>G8II-III</td><td>4780(230)</td><td>—</td><td>2.54(4/9)</td><td>2.5(2)</td></tr><tr><td>φ<sub>2</sub> Ori</td><td>K0IIb</td><td>4670(230)</td><td>2.20(2)</td><td>2.263(6/5)</td><td>2.1(2)</td></tr><tr><td>η Col</td><td>K0III</td><td>4660(230)</td><td>2.48(3)</td><td>2.38(7/8)</td><td>2.7(2)</td></tr><tr><td>λ Hya</td><td>K0III</td><td>4660(230)</td><td>—</td><td>2.62(14/29)</td><td>—</td></tr><tr><td>γ Lib</td><td>K0III</td><td>4660(230)</td><td>—</td><td>2.31(15/15)</td><td>—</td></tr><tr><td>o Sgr</td><td>K0III</td><td>4660(230)</td><td>—</td><td>2.50(12/27)</td><td>2.4(2)</td></tr><tr><td>ℓ Eri</td><td>K0.5IIb</td><td>4600(220)</td><td>2.18(2)</td><td>2.35(4/37)</td><td>2.5(2)</td></tr><tr><td>θ Psc</td><td>K0.5III</td><td>4580(220)</td><td>2.00(2)</td><td>2.11(2/24)</td><td>2.1(2)</td></tr><tr><td>γ Scl</td><td>K1III</td><td>4510(220)</td><td>2.13(3)</td><td>2.039(14/7)</td><td>2.2(2)</td></tr><tr><td>HR 2113</td><td>K1.5III</td><td>4440(220)</td><td>—</td><td>2.48(2)</td><td>2.4(2)</td></tr><tr><td>ε TrA</td><td>K1.5III</td><td>4440(220)</td><td>2.56(7)</td><td>2.43(2)</td><td>—</td></tr><tr><td>η Cet</td><td>K1.5III</td><td>4440(220)</td><td>3.44(4)</td><td>3.323(28/9)</td><td>3.3(2)</td></tr><tr><td>HR 3282</td><td>K2.5II-III</td><td>4330(210)</td><td>2.54(4)</td><td>2.32(7)</td><td>—</td></tr><tr><td>HR 2411</td><td>K3III</td><td>4260(210)</td><td>1.90(3)</td><td>1.76(15/15)</td><td>—</td></tr><tr><td>51 Hya</td><td>K3III</td><td>4260(210)</td><td>2.28(3)</td><td>2.23(7/19)</td><td>2.4(2)</td></tr></tbody></table>

<figure><img src="image_11.png" /><figcaption>**Figure 8.** Value of the angular diameter of each calibrator of our observing sample found in the CCSL, versus the final value deduced from SPIDAST, fitting the MARCS model spectra on calibrated templates of Cohen et al. (1999) (red dots), and on IRAS-LRS measurements (black dots). Solid lines: ±10 per cent thresholds.</figcaption></figure>

each calibrator, and its associated error. Fitting the stellar models of Kurucz (1979) on the absolutely calibrated spectrophotometric templates of Cohen et al. (1999), the CCSL contains 374 calibrators, with LD angular diameters ranging from 1 to 3 mas. Since our selected calibrators might have no angular diameter estimate (one-third of the calibrators of our data set have no value found in the catalogues), or poorly estimated, we have added to SPIDAST various routines for determination of the angular diameter from indirect methods, presented in detail in Cruzalèbes et al. (2010). These routines compute the angular diameter:

(i) by combining the linear diameter, deduced from absolute luminosity and effective temperature from the Morgan Keenan (MK) spectral type, with the parallax;

(ii) by using experimental laws based on the interstellar-corrected colour index (surface brightness method);

(iii) by scaling synthetic spectra on broad-band photometric measurements (infrared flux method); and

(iv) by fitting synthetic spectra on infrared spectrophotometric measurements (see Section 6).

Table 2 gives the angular diameter values of our calibrators found in the CCSL, as well as our final estimates, derived from the fit of MARCS + TURBOSPECTRUM synthetic spectra on photometric or spectrophotometric measurements. In the last column, we also give the estimates found in the JMMC Stellar Diameter Catalogue (JSDC) of Lafrasse et al. (2010), with accuracies between 7 and 10 per cent (Delfosse 2004; Bonneau et al. 2006).

Fig. 8 shows the comparison of our estimates with the CCSL values. The mean difference between our estimates and the CCSL is 5 per cent. The mean difference for the three targets with the available calibrated templates of Cohen (i.e. $\phi_2$ Ori, $\gamma$ Scl and $\eta$ Cet) is lower: 3.5 per cent. This shows the satisfactory agreement between our results and the CCSL, which confirms the reliability of our approach for determination of the calibrator angular diameters, implemented in SPIDAST.

### 4.2.2 Modelling the calibrator visibility and triple product

Once we determined the angular diameter of each calibrator, and its associated uncertainty, we derive the IRF for each SPI observable, computing the ratio of the observable measured on the calibrator, to the synthetic observable, assumed to represent the ‘true’ calibrator observable. Each synthetic SPI observable is derived from the synthetic coherent flux $C_{\lambda}^{\text{model}}$:

(i) the synthetic flux (equation 1) is given by the synthetic coherent flux computed for the null baseline (visibility unity);

(ii) the squared synthetic visibility (equation 2) by the ratio of the squared synthetic coherent flux to the synthetic squared flux;

(iii) the synthetic bispectrum (equation 3) by the product of the synthetic coherent flux for the three baselines composing each triplet; and

(iv) and the synthetic triple product (equation 4) by the ratio of the synthetic bispectrum to the cubed synthetic flux.

Applying the van Cittert-Zernike theorem, the model coherent flux of each calibrator, assumed to emit a luminous intensity with circular symmetry, is derived from the Hankel transform of the source intrinsic spectral radiance $L_{\lambda}$ (in W m$^{-2}$ µm$^{-1}$ sr$^{-1}$)

$$C_{\lambda}^{\text{model}} = \pi \frac{\phi^2}{4} \int_{0}^{\infty} L_{\lambda}(r) J_{0}\left(\pi r \phi \frac{B_{ij}}{\lambda}\right) r dr, \quad (7)$$

where $\phi$ is the stellar angular diameter, $r$ the impact parameter, i.e. the distance from the centre of the stellar disc ($r=0$ at disc centre), $B_{ij}$ the length of the projected baseline of the $(i,j)$ pair of apertures and $J_0$ the Bessel function of the order of zero.

Since compact photospheres are known to deviate from simple uniform discs (see Hajian et al. 1998, and references therein), radially decreasing from their photometric centre, we use the MARCS + TURBOSPECTRUM codes to produce models of Centre-to-Limb Variation (CLV) profiles, with input parameters derived from the MK spectral type. Fig. 9 shows the $\lambda-\mu$ map of typical CLV profiles given by the model, where the reduced radial coordinate $\mu$ is given by $\mu = \sqrt{1-r^2}$ (Young 2003).

With ‘partially’ resolved calibrators, the deviation from the uniform-disc (UD) model increases with the baseline length. Fig. 10
<figure><img src="image_12.png" /><figcaption>**Figure 9.** λ-μ map of CLV profiles produced by the MARCS + TURBOSPEC-TRUM code, in the K band, for a K0III-type star.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>**Figure 10.** Theoretical interferometric observables of the calibrator η Ceti, produced by the LD-MARCS model (in colour), and by the UD model with the same value of the angular diameter (in black). Top panel: squared visibilities for the baselines G1–I1 (blue), D0–G1 (green) and D0–I1 (red). Bottom panel: real part of the triple product produced by the models, for the baseline triplet D0–G1–I1. For each curve, the thickness gives the uncertainty.</figcaption></figure>

shows the squared visibility and the real part of the triple product produced by the MARCS and the UD models, for the calibrator η Ceti ($\phi_{final} = 3.32 \pm 0.02$ mas), deduced from measurements obtained with the D0–I1–G1 baseline triplet (G1–I1 = 46 m, D0–G1 = 69 m and D0–I1 = 79 m). In this example, the differences in visibility and triple product between the UD model and the LD-MARCS model are larger than the uncertainty of the synthetic observables (derived

by propagating the uncertainty on the angular diameter), hence, it is justified to use the LD model instead of the UD model.

### 4.2.3 Correcting for the degradations of the IRF

The main source of non-stationary instabilities which affects the fringe formation is the atmospheric-phase turbulence (Roddier 1981). This issue is partially solved thanks to the use of fringe-tracking servo-loops scanning the fringes with cycle periods smaller than the seeing coherence time, defined as the time over which the phase fluctuations remain coherent (Davis & Tango 1996; Kellerer & Tokovinin 2007). Thus, the residual instabilities are expected to cause only slow drifts of the IRF between observations of the science targets and their associated calibrators. In order to estimate the response at the time when each science target is measured, linear interpolation between the successive measurements on its calibrator(s) is legitimate (Perrin 2003), provided that the ‘calibrator-science-calibrator’ observation sequence repeated several times and rapidly samples the evaluation of the IRF (Berger, Dumas & Kaüfer 2011).

To perform the calibration of the squared visibility and triple product data, SPIDAST applies the following procedure, working at spectral channel level.

(i) Compute the measured IRF for each exposure, given by the ratio of the calibrator measured observable to the calibrator synthetic observable (at the same spectral resolution, for the same baseline or baseline triplet). Values and uncertainties are derived from the fourth-order approximation formulae of Winzer (2000). For triple product data, the approximation formulae apply separately to the real and imaginary parts of the complex ratio. Since the calibrator angular diameter uncertainties are smaller than 10 per cent (see Table 2), the uncertainty of a model observable is simply obtained from its first partial derivative w.r.t. the angular diameter, multiplied by the angular diameter uncertainty.

(ii) Average the IRF measurements over each OB. As weight is associated with each exposure, we use the ratio of the associated weight (defined in Section 3.2) to the variance of the calibrator ‘raw’ observable.

(iii) Determine the IRF at the mean time of each exposure obtained with the science target from linear interpolation, or polynomial fit, on the averaged IRF measurements. Fig. 11 shows the temporal variation of the measured and interpolated response in squared visibility and triple product, obtained during one observing night, with three successive calibrators.

(iv) Divide the observable measured on the science target for each exposure of a given OB, by each interpolated value of the IRF obtained for the same OB (using Winzer’s formulae). Applying a method similar to that used with the *getCal* tool of the Cal. Inst. Tech. (2008), a weight is associated with each calibrated ratio, which combines information on the observations of the calibrator and of the science target (science–calibrator angular separation and observation-time delay, coherence time and number of non-aberrant data).

(v) Compute the weighted average over each OB of the ratios used in the calibration, which gives the final calibrated SPI observable (see Section 4.2.4). The argument of the final calibrated triple product gives the final calibrated closure phase.

Although the calibration procedure used by SPIDAST is based on the standard calibration method, great care has been taken to provide reliable uncertainty estimates, thanks to the use of weights tracing
<figure><img src="image_14.png" /><figcaption>**Figure 11.** Temporal variation of the IRF ($\lambda$ = 2.2 µm) during one observing night. Blue dots: response values averaged over each OB, successively obtained with the calibrators $\gamma$ Lib, $\varepsilon$ TrA and $\phi$ Sgr. Red dots: response values interpolated at the time of observation of the science targets (red dashed lines: interpolated response). From top to bottom: response in squared visibility for the baselines A0–D0, D0–H0, A0–H0, real and imaginary parts of the response in triple product for the triplet A0–D0–H0.</figcaption></figure>

the data quality at different levels of the processing. Thus, SPIDAST produces large uncertainties on final calibrated observables obtained with input data of poor quality. Estimating reliable error bars is of crucial importance to get calibrated data, usable for the fitting process of parametric models, even on data of uneven quality (see Section 6).

### 4.2.4 Comparing the final calibrated observables with AMDLIB

Although the calibration routine provided with the AMDLIB package is not yet officially validated (Malbet et al. 2011), we compute the differences in squared visibility and closure phase between SPIDAST and AMDLIB, and compare them with the associated uncertainties.

Fig. 12 shows the results obtained on the science target $\delta$ Oph, with the baseline triplet A0–D0–H0, under good and poor seeing conditions. The histograms of the uncertainties show that the uncertainties in squared visibility estimated with AMDLIB are smaller than with SPIDAST (whatever the seeing), while it is the contrary for the uncertainties in closure phase.

To explain the disagreement between the uncertainties in squared visibility, we remind that SPIDAST takes into account instrumental and environmental conditions in the calculation of the uncertainties, while AMDLIB uses only unweighted variance. Thus, we consider as underestimated the uncertainties in squared visibility produced by AMDLIB.

Regarding the closure phase, it is more difficult to explain why the uncertainties from AMDLIB are larger than the ones from SPIDAST. Indeed, no details are found in Malbet et al. (2011), expected to describe the AMDLIB calculation. SPIDAST computes the uncertainties on the real and imaginary parts of the calibrated triple product, and derives the uncertainties on the final closure phase from them. A deep analysis based on the source code of the calibration routine provided by AMDLIB reveals that AMDLIB computes the uncertainties in the final closure phase, adding the uncertainties in closure phase of the science target and of the estimated transfer function. Since the closure phase is less stable than the triple product, the final uncertainties in closure phase handled with AMDLIB are larger than with SPIDAST.

To evaluate the agreement between the two procedures, we compare the differences between the final calibrated SPI observables they produce to the associated uncertainties. Since the median differences are smaller than the median uncertainties (Fig. 12), we conclude that the agreement between the two software packages is satisfactory, SPIDAST giving more reliable uncertainties than AMDLIB.

## 5 MEASURING THE DEVIATION FROM CIRCULAR SYMMETRY

The closure phase $\psi_{123}$ couples the phases $\varphi$ of the Fourier transform of the source brightness distribution at the three spatial frequencies $f_1, f_2$ and $f_3 = f_1 + f_2$, probed by the three baselines in the following way: $\psi_{123} = \varphi(f_1) + \varphi(f_2) - \varphi(f_3)$. In the case of a centrosymmetrical brightness distribution on-axis, the phase of the Fourier transform is uniformly zero and $\psi_{123}$ is naturally equal to zero. If this source is off-axis, the phase is linear w.r.t. the spatial frequency, so that $\varphi(f_3) = \varphi(f_1) + \varphi(f_2)$, and, subsequently, also here $\psi_{123}$ is equal to zero. Thus, a non-zero phase closure is a clear indication, at least qualitatively, for a deviation from circular symmetry.

In Section 5.2, we introduce a new parameter (called the centrosymmetry parameter, CSP), based on the triple product, more sensitive to deviation from centrosymmetry than the closure phase.

### 5.1 Global closure phase

The integration of the real and imaginary parts of the triple product $\mathcal{T}_{\text{true}}$, over the observation spectral band [$\lambda_{\text{min}}; \lambda_{\text{max}}$], leads to the global closure phase, defined by Ragland et al. (2006) and Tatebe et al. (2006) as

$$ \psi_{\text{band}} = \operatorname{atan} \frac{\int_{\lambda_{\text{min}}}^{\lambda_{\text{max}}} \mathfrak{T}_{\text{true}}(\lambda) \, d\lambda}{\int_{\lambda_{\text{min}}}^{\lambda_{\text{max}}} \mathfrak{T}_{\text{true}}(\lambda) \, d\lambda}. \qquad (8) $$

A value of the global closure phase close to zero (less than 1 $^{\circ}$, for our sample) suggests a high degree of centrosymmetry of the brightness distribution, related to the observation spectral band. We compute the uncertainty on the global closure phase thanks to the direct bootstrap method, using random sampling with replacement of the spectral data set of the calibrated triple product. The values of the global closure phases in K, and their associated uncertainties, are given in the third column of Table 3, for each science target.

### 5.2 Centrosymmetry parameter

When the positive and negative values of $\mathfrak{T}_{\text{true}}(\lambda)$ along the spectral domain almost mutually compensate (as seen on Fig. 13 for TX
<figure><img src="image_15.png" /><figcaption>**Figure 12.** Histograms of calibrated quantities pertaining to the science target $\delta$ Oph, obtained in good (blue) and poor (red) seeing conditions. From top to bottom: squared visibilities with the baselines A0–D0, D0–H0, A0–H0 and closure phase with the triplet A0–D0–H0. Long panels of the first column: difference between SPIDAST and AMDLIB. Small panels of the other columns: associated uncertainties with good (blue: SPIDAST, turquoise: AMDLIB) and poor seeing (red: SPIDAST, pink: AMDLIB).</figcaption></figure>

<table><caption>**Table 3.** CSP values and global closure phases for the science targets. The last column gives the ratio of the associated SNRs.</caption><thead><tr><th>Name</th><th>CSP (°)</th><th>ψ<sub>band</sub> (°)</th><th>SNR<sub>CSP</sub>/SNR<sub>ψ</sub></th></tr></thead><tbody><tr><td>δ Oph</td><td>0.76(3)</td><td>0.16(4)</td><td>7.0</td></tr><tr><td>α Car</td><td>0.84(3)</td><td>-0.22(5)</td><td>5.7</td></tr><tr><td>L<sub>2</sub> Pup</td><td>0.90(4)</td><td>0.39(5)</td><td>3.2</td></tr><tr><td>β Cet</td><td>0.98(4)</td><td>0.32(6)</td><td>4.6</td></tr><tr><td>ζ Ara</td><td>1.12(4)</td><td>-0.22(6)</td><td>8.1</td></tr><tr><td>α TrA</td><td>1.21(5)</td><td>0.20(7)</td><td>8.3</td></tr><tr><td>α Hya</td><td>1.21(9)</td><td>-0.16(11)</td><td>9.4</td></tr><tr><td>TW Oph</td><td>2.94(5)</td><td>-2.91(5)</td><td>1.1</td></tr><tr><td>CE Tau</td><td>3.01(11)</td><td>-0.85(17)</td><td>5.5</td></tr><tr><td>γ Hyi</td><td>3.25(11)</td><td>2.83(13)</td><td>1.4</td></tr><tr><td>o<sub>1</sub> CMa</td><td>4.41(197)</td><td>2.88(155)</td><td>1.2</td></tr><tr><td>σ Lib</td><td>5.12(20)</td><td>-1.90(33)</td><td>4.4</td></tr><tr><td>γ Ret</td><td>5.16(37)</td><td>-4.44(115)</td><td>3.6</td></tr><tr><td>TX Psc</td><td>5.34(25)</td><td>-2.17(32)</td><td>3.4</td></tr><tr><td>o<sub>1</sub> Ori</td><td>8.10(82)</td><td>7.48(84)</td><td>1.1</td></tr><tr><td>W Ori</td><td>10.58(46)</td><td>1.35(79)</td><td>13.3</td></tr><tr><td>R Scl</td><td>11.63(59)</td><td>-10.63(67)</td><td>1.2</td></tr><tr><td>T Cet</td><td>27.88(51)</td><td>-28.79(41)</td><td>0.8</td></tr></tbody></table>

Psc and W Ori), the spectral integration produces a nearly null global closure phase. As mentioned above, this can be taken as a hint for centrosymmetry. However, $\mathfrak{T}_{true}(\lambda)$ is clearly non-null in some parts of the spectrum, and this is a hint for a deviation from

centrosymmetry. To rule out the contradiction, and so as to consider this latter possibility, we introduce a new estimator, that we call the CSP, similar to the global closure phase, but using instead the absolute value of $\mathfrak{T}_{true}(\lambda)$ in the numerator, and the modulus of $\mathfrak{T}_{true}(\lambda)$ in the denominator,

$$ \text{CSP} = \text{asin} \frac{\int_{\lambda_{\text{min}}}^{\lambda_{\text{max}}} |\mathfrak{T}_{\text{true}}(\lambda)| \, d\lambda}{\int_{\lambda_{\text{min}}}^{\lambda_{\text{max}}} |\mathfrak{T}_{\text{true}}(\lambda)| \, d\lambda}. \quad (9) $$

As for the global closure phase, a small CSP value (less than 2°, for our sample) suggests a high degree of centrosymmetry. Significantly high CSP values, as judged from their uncertainties (bootstrap method again), require to use asymmetric brightness distribution models, for the fitting on the SPI data. The values of the global CSP in the K band, and their uncertainties, are given in the second column of Table 3.

To illustrate the difference between the global closure phase and the CSP, we plot in Fig. 13 three quantities, for three scientific targets showing a global closure phase close to zero,

(i) in the panels on the left is shown the imaginary part of the triple product $\mathfrak{T}_{true}$, divided by the spectral mean<fnref n="4" /> of the real part

<footnote n="4">⁴ Defined as the integral w.r.t. the wavelength, divided by the spectral bandwidth $\Delta\lambda = \lambda_{max} - \lambda_{min}$.</footnote>
<figure><img src="image_16.png" /><figcaption>**Figure 13.** Left-hand panels: spectral variation of the imaginary part of the triple product, divided by the spectral average of the real part, for three science targets (the 'spectral mean' quoted above each panel corresponds to the spectral average of the displayed quantity; for the left-hand panels, it thus equals the tangent of the global closure phase). Central panels: absolute value of the imaginary part of the triple product, divided by the spectral average of the modulus (for the central panel, the 'spectral mean' equals the sinus of the CSP). Right-hand panels: absolute value of the real part of the triple product, divided by the spectral average of the modulus.</figcaption></figure>

⟨ϑT<sub>true</sub>⟩. Note that the spectral mean of this ratio is equal to the tangent of the global closure phase, as defined in equation (8);

(ii) in the panels forming the central column is shown the absolute value of the imaginary part of the triple product |ϑT<sub>true</sub>|, divided by the spectral mean of the modulus ⟨|T<sub>true</sub>|⟩. Note that the spectral mean of this ratio is equal to the sine of the global CSP, as defined in equation (9); and

(iii) in the panels on the right is shown the absolute value of the real part of the triple product |ϑT<sub>true</sub>|, divided by the spectral mean of the modulus ⟨|T<sub>true</sub>|⟩.

The three panels on the left show three typical behaviours for
ϑT<sub>true</sub> :

(i) uniformly close to zero (δ Oph);

(ii) decreasing symmetrically around zero (TX Psc); and

(iii) increasing symmetrically around zero (W Ori).

For each of these situations, the integration over the spectrum pro-
duces a nearly null result (hint for centrosymmetry). However, the
CSP allows to detect deviations from centrosymmetry (see Fig. 14).

If the real and imaginary parts of the triple product are wavelength
independent, for the observation spectral range, one can show that
sin CSP = |ψ<sub>band</sub>|. If not, no analytical relation linking the two
quantities can be derived from equations (8) and (9), because of the
integrals in the numerators and the denominators. Fig. 14 shows
the CSP in *K* displayed versus the global closure phase in the same
spectral band, for each science target, choosing the OB which gives

the smallest uncertainty on the CSP. Almost all the stars fall along
the diagonals of the (ψ<sub>band</sub>, CSP) diagram, which trace the relation
CSP = |ψ<sub>band</sub>|. One (W Ori), however, is flagged as asymmetrical
with the CSP indicator, but not with ψ<sub>band</sub>. This is precisely why
we favour the CSP over ψ<sub>band</sub>, as discussed above. T Cet has the
largest values for CSP and ψ<sub>band</sub> as well, indicating large deviation
from centrosymmetry (may be due to strong asymmetries at the
wavelengths of the CO bands or even suggestive of the presence of
a binary companion). The second largest CSP value pertains to R
Scl, recently revealed as a wide binary by Maercker et al. (2012),
using the ALMA array at millimetric wavelengths. In addition, the
last column of Table 3 shows that, except for T Cet, the SNR for the
CSP is comparable to or larger than the SNR for the global closure
phase.

6 FITTING PARAMETRIC CHROMATIC MODELS

In order to interpret the final calibrated SPI observables, SPIDAST
provides a fitting routine, based on the modified gradient-expansion
algorithm, very similar to the algorithm invented by Levenberg
(1944), and improved by Marquardt (1963). The fit applies on any
calibrated SPI measurements (to be chosen between visibility, flux,
coherent flux, bispectrum, triple product and closure phase), using
a library of single-component or composite parametric chromatic
models, characterized by the Fourier spectrum of their intensity dis-
tribution, and the associated first-order partial derivatives w.r.t. the
<figure><img src="image_17.png" /><figcaption>**Figure 14.** The CSP versus $\psi_K$, the global closure phase in $K$, for the science targets. The red inset in the upper-right corner enlarges the group of points located at the bottom, around the null value of the global closure phase, pertaining to the targets $\delta$ Oph, $\alpha$ Car, $L_2$ Pup, $\beta$ Cet, $\zeta$ Ara, $\alpha$ TrA and $\alpha$ Hya. The dash-dotted diagonals trace the relation CSP = $|\psi_K|$.</figcaption></figure>

model parameters. Fig. 15 shows examples of visibility-fitting results obtained with SPIDAST, using two models (UD and LD-MARCS). This fitting routine, applied on spectrophotometric data, is used to determine the calibrator angular diameters from the MARCS + TURBOSPECTRUM synthetic spectra (see Section 4.2.1).

Note that the LITPRO software,<fnref n="5" /> developed by the JMMC working group ‘Model fitting’, uses a set of elementary geometrical and centre-to-limb darkening functions as well (Tallon-Bosc et al. 2008). However, it does not offer, contrary to SPIDAST, the possibility to fit stellar-disc models with synthetic tabulated radiance data (or exitance, for fits on spectrophotometric measurements), such as those produced by the MARCS + TURBOSPECTRUM code.

Since the uncertainties of the final calibrated data are not normally distributed, the covariance matrix, that comes out of the $\chi^2$ fit, cannot be used to infer the parameter uncertainties (Press et al. 2007). To compute the uncertainties, SPIDAST uses the residual-bootstrap method, described in detail in Cruzalèbes et al. (2010):

(i) ‘synthetic’ data sets are produced from random sampling with replacement of the Pearson residuals (difference between calibrated and fitted values divided by the uncertainty of the observed value), added to the initial fitted values;

(ii) fitting the model on these new data sets produces a set of $\chi^2$ minima, following a probability distribution, from which we extract the boundaries of the confidence interval, with a given confidence level; and

(iii) the parameter values associated with these $\chi^2$ boundaries give the upper and lower limits of the parameter estimates, leading to asymmetric uncertainties.

The whole procedure allows the computing of reliable estimates of the parameter uncertainties.

## 7 CONCLUSION

In this paper, we introduce the new SPIDAST software, developed since 2006 with the aim to reduce, calibrate and interpret the visibility and triple product measurements obtained with the VLTI/AMBER facility. SPIDAST contains a whole set of modules, which can be launched separately or in an automatic batch file, and summarized hereafter.

(i) The raw data reduction used by SPIDAST computes the weighted average of non-aberrant data, at spectral-channel level, providing estimates of the SPI observables using an automated procedure, while the method presently used with AMDLIB selects the ‘best’ frames according to a quality threshold determined a posteriori, after several trials.

(ii) The wavelength calibration procedure performed by SPIDAST provides spectral shifts following a polynomial law, tracing them at the channel level. This is done by computing the cross-correlation of the measured spectra of the calibrators with their synthetic spectra produced by the MARCS model, while AMDLIB only provides a constant spectral shift over the spectrum, from the correlation with the telluric lines.

(iii) For selected calibrators not included in the calibrator catalogues, SPIDAST provides several routines for estimation of the angular diameter with indirect methods, the most accurate being the fit of stellar-atmosphere model spectra given by MARCS on spectrophotometric data. The calibrator synthetic observables are derived using CLV functions produced by MARCS.

(iv) To obtain an accurate interferometric calibration (via an automated procedure), SPIDAST: (1) divides the calibrator raw data with the associated synthetic observables in each spectral channel, which gives the instrumental response function in squared visibility and triple product; (2) interpolates or fits the response at the time of each exposure on the science target; (3) divides the science raw data with the interpolated/fitted response, which gives the science calibrated observables for each exposure; and (4) computes their weighted average over each OB. At each processing step, the uncertainties are computed, thanks to the bootstrap method applied on the weighted means.

(v) Using the real and imaginary parts of the calibrated triple product, SPIDAST measures the deviation from centrosymmetry of the brightness distribution of each science target in the observation spectral band, thanks to a new parameter, more sensitive to asymmetries than the global closure phase.

(vi) Finally, SPIDAST proposes a complete fitting tool, using a set of parametric and chromatic models, and accepting input tables of flux/intensity synthetic data.

Such a careful calibration process of SPI data is a crucial step for their trustworthy astrophysical exploitation, which is the topic of two associated papers (Cruzalèbes et al., in preparation). Parameter extraction using non-linear fits of source model, as well as aperture-synthesis image reconstruction, need reliable estimates of the calibrated observables, with robust uncertainties. Our reduction, calibration and fitting routines also apply to any other spectral data sets, including spectroscopic data. We made the SPIDAST software public:<fnref n="6" /> the source code of any program of our software suite can be obtained by sending an e-mail to the first author of this paper.

<footnote n="5">⁵ www.jmmc.fr/litpro_page.htm</footnote>

<footnote n="6">⁶ https://forge.oca.eu/trac/spidast</footnote>
<figure><img src="image_18.png" /><figcaption>**Figure 15.** Fits on the visibility data of δ Oph obtained with the baselines A0–D0, D0–H0 and A0–H0, for one single OB. Left-hand panel: fit of the UD model; right-hand panel: fit of the LD-MARCS model (models in dashed lines). Calibrated data are shown without error bars, for clarity purpose. φ is the best-fitting angular outer diameter (in mas), χ²_R the reduced χ² and R2 the adjusted coefficient of determination. Bottom panels: residuals (calibrated data-fitted values).</figcaption></figure>

## ACKNOWLEDGEMENTS

The authors thank the ESO-Paranal VLTI team for supporting their AMBER observations, especially the night astronomers A. Mérand, G. Montagnier, F. Patru, J.-B. Le Bouquin, S. Rengaswamy, and W. J. de Wit, the VLTI group coordinator S. Brillant, and the telescope and instrument operators A. and J. Cortes, A. Pino, C. Herrera, D. Castex, S. Cerda, and C. Cid. They also thank the Programme National de Physique Stellaire (PNPS) for supporting part of this collaborative research. AJ is grateful to B. Plez, K. Eriksson and T. Masseron for their ongoing support on the use of the MARCS code. SS was (partly) supported by the Austrian Science Fund through FWF project P19503-N16; AC is post-doctoral fellow from F.R.S.-FNRS (Belgium; grant 2.4513.11); and EP is supported by PRODEX. The data were partially reduced using the publicly available data reduction software package AMDLIB, kindly provided by the Jean-Marie Mariotti Center. This study used the SIMBAD and VIZIER data bases at the CDS, Strasbourg (France), and NASA ADS bibliographic services.

## REFERENCES

* Alvarez R., Plez B., 1998, A&A, 330, 1109

* Berger J.-P., Dumas C., Kaüfer A., 2011, VLT Paranal Science Operations AMBER User Manual. ESO, Garching

* Boden A. F., 2007, New Astron., 51, 617

* Bonneau D. et al., 2006, A&A, 456, 789

* Bordé P., Coudé du Foresto V., Chagnon G., Perrin G., 2002, A&A, 393, 183

* Cal. Inst. Tech. 2008, The GetCal Interferometric Observation Planning Tool Suite, 2.10.7 edn. NASA Exoplanet Science Inst., Pasadena, CA

* Cohen M., Walker R. G., Carter B., Hammersley P., Kidger M., Noguchi K., 1999, AJ, 117, 1864

* Cohen M., Wheaton W. A., Megeath S. T., 2003, AJ, 126, 1090

* Cruzalèbes P., Spang A., Sacuto S., 2008, in Kaufer A., Kerber F., eds, The 2007 ESO Instrument Calibration Workshop, ESO Astrophysics Symposia European Southern Observatory, Calibration of AMBER Visibilities at Low Spectral Resolution. Springer-Verlag, Berlin, p. 479

* Cruzalèbes P., Jorissen A., Sacuto S., Bonneau D., 2010, A&A, 515, A6

* Davis J., Tango W., 1996, PASP, 108, 456

* Decin L., Eriksson K., 2007, A&A, 472, 1041

* Decin L., Morris P. W., Appleton P. N., Charmandaris V., Armus L., Houck J. R., 2004, ApJS, 154, 408

* Delfosse X., 2004, Technical report, Stellar Diameter Estimation from Photospheric Indices. JMMC, Grenoble, France

* Domiciano de Souza A., Bendjoya P., Vakili F., Millour F., Petrov R. G., 2008, A&A, 489, L5

* Efron B., 1979, Ann. Stat., 7, 1

* Efron B., Tibshirani R., 1993, Chapman & Hall/CRC Monographs on Statistics & Applied Probability: An Introduction to the Bootstrap. Chapman & Hall/CRC Press, Boca Raton

* Guizar-Sicairos M., Thurman S. T., Fienup J. R., 2008, Opt. Lett., 33, 156

* Gustafsson B., Edvardsson B., Eriksson K., Jørgensen U. G., Nordlund Å., Plez B., 2008, A&A, 486, 951

* Hajian A. R. et al., 1998, ApJ, 496, 484

* Hoaglin D. C., Mosteller F., Tukey J. W., 1983, Understanding Robust and Exploratory Data Analysis. Wiley, New York

* Kellerer A., Tokovinin A., 2007, A&A, 461, 775

* Kraus S. et al., 2009, A&A, 497, 195

* Kurucz R. L., 1979, ApJS, 40, 1

* Lafrasse S., Mella G., Bonneau D., Duvert G., Delfosse X., Chesneau O., Chelli A., 2010, in Danchi W. C., Delplancke F., Rajagopal J. K., eds, SPIE Conf. Ser. Vol. 7734, Optical and Infrared Interferometry II. SPIE, Bellingham, p. 77344E

* Levenberg K., 1944, Q. Appl. Math., 2, 164

* Maercker M. et al., 2012, Nat, 490, 232

* Malbet F., Millour F., Absil O., Duvert G., 2011, User Manual – Issue 3.0.3 JMMC-MAN-2720-0001, AMBER Data Reduction Software, JMMC, available at www.jmmc.fr/doc/approved/JMMC-MAN-2720-0001.pdf

* Marquardt D. W., 1963, Indian J. Indust. Appl. Math., 11, 431

* Martí-Vidal I., Marcaide J. M., Quirrenbach A., Ohnaka K., Guirado J. C., Wittkowski M., 2011, A&A, 529, A115

* Mérand A., Bordé P., Coudé du Foresto V., 2005, A&A, 433, 1155

* Mérand A., Stefl S., Bourget P., Ramirez A., Patru F., Haguenauer P., Brillant S., 2010, in Danchi W. C., Delplancke F., Rajagopal J. K., eds, SPIE Conf. Ser. Vol. 7734, Optical and Infrared Interferometry II. SPIE, Bellingham, p. 77342B

* Millour F. et al., 2007, A&A, 464, 107

* Millour F. et al., 2008, in Kaufer A., Kerber F., eds, The 2007 ESO Instrument Calibration Workshop, ESO Astrophysics Symposia European Southern Observatory, AMBER on the VLTI Data Processing and Calibration Issues. Springer-Verlag, Berlin, p. 461

* Ohnaka K. et al., 2011, A&A, 529, A163

* Perrin G., 2003, A&A, 400, 1173

* Petrov R. G. et al., 2007, A&A, 464, 1
Plez B., 2003, in Munari U., ed., ASP Conf. Ser. Vol. 298, Cool Star Atmospheres and Spectra for GAIA: MARCS Models. Astron. Soc. Pac., San Francisco, p. 189

Plez B., 2012, Astrophys. Source Code Libr., record ascl:1205.004

Press W. H., Teukolsky S. A., Vetterling W. T., Flannery B. P., 2007, Numerical Recipes: The Art of Scientific Computing, 3rd edn. Cambridge Univ. Press, Cambridge

Ragland S. et al., 2006, ApJ, 652, 650

Rix S. A., Pettini M., Leitherer C., Bresolin F., Kudritzki R.-P., Steidel C. C., 2004, ApJ, 615, 98

Robbe-Dubois S. et al., 2007, A&A, 464, 13

Roddier F., 1981, Prog. Opt., 19, 281

Stefl S., Mérand A., LeBouquin J.-B., Rantakyro F., Bourget P., Carciofin A., Rivinius T., Baade D., 2011, in The Origin and Fate of the Sun: Evolution of Solar-mass Stars Observed with High Angular Resolution Problems of the AMBER Wavelength Calibration. ESO, Garching

Tallon-Bosc I. et al., 2008, in Schöller M., Danchi W. C., Delplancke F., eds., SPIE Conf. Ser. Vol. 7013, Optical and Infrared Interferometry. SPIE, Belligham, p. 70131J

Tatebe K., Chandler A. A., Hale D. D. S., Townes C. H., 2006, ApJ, 652, 666

Tatulli E. et al., 2007, A&A, 464, 29

van Belle G. T., van Belle G., 2005, PASP, 117, 1263

Vera E., Torres S., 2008, in EURASIP, ed., 16th ESP Conf. (EUSIPCO 2008) Subpixel Accuracy Analysis of Phase Correlation Registration Methods using Aliased Imagery. EURASIP, Lausanne, Switzerland

Weigelt G. P., 1977, Opt. Commun., 21, 55

Weigelt G. et al., 2011, A&A, 527, A103

Winzer P., 2000, Rev. Sci. Instrum., 71, 1447

Wittkowski M., Boboltz D. A., Driebe T., Le Bouquin J.-B., Millour F., Ohnaka K., Scholz M., 2008, A&A, 479, L21

Young J., 2003, in Perrin G., Malbet F., eds, EAS Publ. Ser. Vol. 6, From Visibilities to Science with Simple Models. EDP Sciences, Les Ulis, France, p. 181

This paper has been typeset from a TeX/LATEX file prepared by the author.