Issue 
A&A
Volume 525, January 2011



Article Number  A21  
Number of page(s)  16  
Section  Numerical methods and codes  
DOI  https://doi.org/10.1051/00046361/201014511  
Published online  29 November 2010 
Fast raytracing algorithm for circumstellar structures (FRACS)
I. Algorithm description and parameterspace study for midIR interferometry of B[e] stars
UMR 6525 CNRS, H. FIZEAUUniversité de Nice SophiaAntipolis, Observatoire
de la Côte d’Azur,
Campus Valrose,
06108
Nice Cedex 2,
France
email: nicolin@unice.fr
Received:
26
March
2010
Accepted:
9
September
2010
Aims. The physical interpretation of spectrointerferometric data is strongly modeldependent. On one hand, models involving elaborate radiative transfer solvers are too time consuming in general to perform an automatic fitting procedure and derive astrophysical quantities and their related errors. On the other hand, using simple geometrical models does not give sufficient insights into the physics of the object. We propose to stand in between these two extreme approaches by using a physical but still simple parameterised model for the object under consideration. Based on this philosophy, we developed a numerical tool optimised for midinfrared (midIR) interferometry, the fast raytracing algorithm for circumstellar structures (FRACS), which can be used as a standalone model, or as an aid for a more advanced physical description or even for elaborating observation strategies.
Methods. FRACS is based on the raytracing technique without scattering, but supplemented with the use of quadtree meshes and the full symmetries of the axisymmetrical problem to significantly decrease the necessary computing time to obtain e.g. monochromatic images and visibilities. We applied FRACS in a theoretical study of the dusty circumstellar environments (CSEs) of B[e] supergiants (sgB[e]) in order to determine which information (physical parameters) can be retrieved from present midIR interferometry (flux and visibility).
Results. From a set of selected dusty CSE models typical of sgB[e] stars we show that together with the geometrical parameters (position angle, inclination, inner radius), the temperature structure (inner dust temperature and gradient) can be well constrained by the midIR data alone. Our results also indicate that the determination of the parameters characterising the CSE density structure is more challenging but, in some cases, upper limits as well as correlations on the parameters characterising the mass loss can be obtained. Good constraints for the sgB[e] central continuum emission (central star and inner gas emissions) can be obtained whenever its contribution to the total midIR flux is only as high as a few percents. Raytracing parameterised models such as FRACS are thus well adapted to prepare and/or interpret long wavelengths (from midIR to radio) observations at present (e.g. VLTI/MIDI) and nearfuture (e.g. VLTI/MATISSE, ALMA) interferometers.
Key words: methods: numerical / methods: observational / techniques: high angular resolution / techniques: interferometric / stars: mass loss / stars: emissionline, Be
© ESO, 2010
1. Introduction
When dealing with optical/IR interferometric data, one needs to invoke a model for the understanding of the astrophysical object under consideration. This is because of (1) the low coverage of the uvplane and most of the time because of the lack of the visibility phase; and (2) because our aim is to extract physical parameters from the data. This is particularly true for the MidInfrared Interferometric Instrument (MIDI, Leinert et al. 2003) at the Very Large Telescope Interferometer (VLTI), on which our considerations will be focused. Some pure geometrical information can be recovered through a simple toy model such as Gaussians (see e.g. Leinert et al. 2004; Domiciano de Souza et al. 2007).
However, this approach does not give any insights into the physical nature of the object. One would dream of having a fully consistent model to characterise the object under inspection. In many cases, if not all, a fully consistent model is out of reach and one uses at least a consistent treatment of the radiative transfer. Models based for instance on the Monte Carlo method are very popular (see e.g. Ohnaka et al. 2006; Niccolini & Alcolea 2006; Wolf et al. 1999) for this purpose. Still, the medium density needs to be parameterised and it is not determined in a selfconsistent way. For massive stars for instance, it would be necessary to take into account nonLTE effects including both gas and dust emission of the circumstellar material as well as a full treatment of radiation hydrodynamics. Fitting interferometric data this way is as yet impossible because of computing time limitations.
Of course, solving at least the radiative transfer in a selfconsistent way is already very demanding for the computational resources. Consequently, model parameters cannot be determined in a fully automatic way and the model fitting process must be carried out mostly by hand, or automatised by systematically exploring the parameterspace , the “chibyeye” approach mentioned in Press et al. (1992). The followers of this approach consider the “best fitting” model as their best attempt: a model that is compatible with the data. It is admittedly not perfect, but it is in most of the cases the best that can be done given the difficulty of the task. It is remarkable that a thorough χ^{2} analysis of VLTI/MIDI data of the Herbig Ae star AB Aurigae has been performed by di Folco et al. (2009) which remains to date one of the most achieved studies of this kind. From the χ^{2} analysis, formal errors can be derived and at least the information concerning the constraints for the physical parameters can be quantified. Qualitative information about the correlation of parameters can be pointed out.
The next step after the toy models for the physical characterisation of the astrophysical objects can be made from the pure geometrical model towards the selfconsistency by including and parameterising the object emissivity in the analysis. For instance, Lachaume et al. (2007) and Malbet et al. (2005) use optically thick (i.e. emitting as black bodies) and infinitely thin discs to model the circumstellar environment of premainsequence and B[e] stars. Of course this approach has some restriction when modelling a disc: for instance it cannot handle nearly edgeon disc and an optically thin situation.
We propose an intermediate approach: between the use of simple geometrical models and sophisticated radiative transfer solvers. Indeed, it is a step backwards from the “selfconsistent” radiative transfer treatment, which is in most cases too advanced with regard to the information provided by the interferometric data. For this intermediary approach, we assume a prescribed and parameterised emissivity for the medium. Our purpose is to derive the physical parameters that characterise this emissivity. In the process, we compute intensity maps and most particularly visibility curves from the knowledge of the medium emissivity with a fast raytracing technique (a few seconds depending on map resolution), taking into account the particular symmetries of a disc configuration. Then, the model fitting process can be undertaken in an automatic way with standard methods (see e.g. Levenberg 1944; Marquardt 1963). The techniques we present are designed to be quite general and not tailored to any particular emissivity except for the assumed axisymmetry of the problem under consideration.
Our purpose is twofold. On one hand – as already mentioned – we aim to estimate physical parameters and their errors characterising the circumstellar dusty medium under consideration with as few restrictive assumptions as possible; at least within the obvious limitations of the present model. On the other hand our purpose is to provide the user of a more detailed model, such as a Monte Carlo radiative transfer code, with a first characterisation of the circumstellar matter to start with.
In Sect. 2 we describe the general framework of the proposed ray tracing technique. In particular how to derive the observable from the astrophysical object emissivity. In Sect. 3 we describe the numerical aspects that are specific to the present raytracing technique. In particular, the use of a quadtree mesh and the symmetries that allow us to speed up the computation are detailed. In Sect. 4 we focus our attention on the circumstellar disc of B[e] stars and describe a parametric model of the circumstellar environment. In Sect. 5 we analyse artificial interferometric data generated both from the parametric model itself and from a Monte Carlo radiative transfer code (Niccolini & Alcolea 2006). Our purpose is not to fit any particular object, but to present our guideline to the following question: which physical information can we get from the data? A discussion of our results and the conclusions of our work are given in Sects. 6 and 7 respectively.
2. The raytracing technique
We describe here the FRACS algorithm, developed to study stars with CSEs from midIR interferometric observables (e.g. visibilities, fluxes, closure phases). Although FRACS could be extended to investigate any 3D CSE structures, we focus here on the particular case of axisymmetrical dusty CSEs. This study is motivated by the typical data one can obtain from disclike CSE observed with MIDI, the midIR 2telescope beamcombiner instrument of ESO’s VLTI (Leinert et al. 2003).
2.1. Intensity map
Intensity maps of the object are the primary outputs of the model that we need to compute the visibilities and fluxes that are directly compared to the observations. For this purpose, we integrate the radiation transfer equation along a set of rays (raytracing technique) making use of the symmetries of the problem (see Sect. 3 for details).
The unit vector along the line of sight is given by , i being the inclination between the zaxis and the line of sight and , ŷ the unit vector along the x et yaxis of a cartesian system of coordinates (see Fig. 1), referred to as the “model system” below. The problem is assumed to be invariant by rotation around the zaxis. We define a fictitious image plane by giving two unit vectors Ŷ = −ŷ cosi + ẑ sini and . This particular choice is made making use of the axisymmetry of the problem. Note that for this particular coordinate system (X,Y) the disc position angle (whenever i ≠ 0) is always defined as 90°. The actual image plane, with the Y′ and X′ axis corresponding respectively to North and East, is obtained by rotating the axis of our fictitious image plane by an angle , where PA_{d} is the position angle of the disc with respect to North.
Fig. 1 Coordinate systems. The shaded ellipse represents a disc viewed by the observer. 

Open with DEXTER 
We neglect the scattering of the radiation by dust grains, optimising our approach to long wavelengths (from midIR to radio). This assumption simplifies the radiative transfer equation by removing the scattering term.
We obtain the intensity map at position (X,Y) in the image plane (inclined by i) and at wavelength λ by integrating the transfer equation along the particular ray that passes through the considered point of the image plane. Defining r_{s}(X,Y,i) (simply r_{s} for short) as the position vector along a ray, given in the model system of coordinates by (2)and by introducing the optical depth at wavelength λ and position s along the ray by (3)we obtain (4)where the extinction coefficient because scattering is neglected.
We assume that the CSE is confined within a sphere of radius R_{out}, s varies consequently from to (R^{2} = X^{2} + Y^{2}) in Eq. (4) and in the definition of a ray Eq. (2). This hypothesis can be relaxed without altering the present considerations and the domain of integration of Eq. (4) suitably chosen.
If some radiation sources (e.g. black body spheres) are included in the analysis, an additional term must be added in Eq. (4) whenever a particular ray intersect a source. For a source with specific intensity this additional term is given by , s^{(s)} being the distance at which the ray given by X, Y and i (see Eq. (2)) intersects the outermost (along the ray) source boundary. In that case the lower integration limit in Eq. (4), that is , must also be replaced by s^{(s)}.
2.2. Interferometric observables
From the monochromatic intensity maps at wavelength λ (Eq. (4)) we obtain both the observed fluxes F_{λ} and visibilities V_{λ} for an object at distance d, (5)and (6)where V_{λ} is obtained for a given baseline specified by its projected length B (on the sky, i.e. (X′,Y′) coordinates) and its polar angle PA from North to East (direction of the Y′ axis). Δ and j represent, respectively, PA_{d} − PA and .
3. Numerical considerations
We seek to produce intensity maps within seconds^{1} and we aim for our numerical method to be sufficiently general in order to deal with a large range of density and temperature structures. Given these two relatively tight constraints, the numerical integration of Eq. (4) is not straightforward.
For example we have tested that the 5th order RungeKutta integrators of Press et al. (1992) with adaptive stepsize (as discussed in Steinacker et al. 2006) doest not suit our constraints. Indeed, the step adaption leads to difficulties if sharp edges (e.g. inner cavities) are present in the medium emissivity.
3.1. Mesh generation
Regarding the above mentioned constraints and the different numerical approaches tested, we found that Eq. (4) is more efficiently computed with an adaptive mesh based on a tree data structure (quadtrees/octrees). The mesh purpose is twofold: first, it must guide the computation of Eq. (4) and distribute the integration points along the rays according to the variations of the medium emissivity; second, within the restriction of axissymmetrical situations, the mesh must handle any kind of emissivity. Quad/octree meshes are extensively used in Monte Carlo radiative transfer codes (e.g. see Bianchi 2008; Niccolini & Alcolea 2006; Jonsson 2006; Wolf et al. 1999); the mesh generation algorithm is thoroughly described in Kurosawa & Hillier (2001).
The mesh we use is a cartesian quadtree. Cartesian refers here to the mesh type and not to the system of coordinates we use. Indeed, the mesh is implemented as a nested squared domain (cells) in the ρ − z plane (). The whole mesh is enclosed by the largest cell (the root cell in the tree hierarchy) of size R_{out} in ρ and z . The underlying physical coordinate system is cylindrical (with z > 0) and the mesh cells correspond to a set of two (for z > 0 and z < 0) tori, which are the actual physical volumes.
The mesh generation algorithm consists in recursively dividing each cell in four child cells until the following conditions are simultaneously fulfilled for each cell in the mesh (see Kurosawa & Hillier 2001, for more details): where V_{ξ} is the volume of cell ξ, V_{tot} is the volume of the root cell and α, β and η are parameters controlling the mesh refinement.
In the present work α and β have been fixed to 1, but higher values can be useful for some particular situations where the generated mesh must be tighter than the mesh generated directly from the and T variations. Typically, these situations show up for high optical depths (in this paper, optical depth values do not exceed ≃1 at 10 μm along the rays). The practical choice of α, β and η is obtained from a compromise between execution speed and numerical accuracy of the Eq. (4) integration (e.g. typical values of η range from 10^{5} to 10^{4}). When dealing with optically thick situations, a supplementary conditions can be added to Eqs. (7) and (8) in order to prescribe an upper limit to the cell optical depth. For instance, making use of the computation of the integral in Eq. (7), one can add the following criterion for cell ξ (whose centre is (ρ_{ξ},z_{ξ}) and size Δ_{ξ}) (9)where Δτ_{lim} is the prescibed upper limit to the cell optical depth.
For the moderate optical depths reached in this work, with values of η down to 10^{5} and Δτ_{lim} set to 10^{2}, the criteria of Eqs. (7) and (8) are the leading conditions to the mesh refinement.
Figure 2 shows the mesh obtained in the particular case of a B[e] circumstellar disc (see Sect. 4) for models whose parameters are given in Table 2 (see caption for more details).
The volume integrals in Eqs. (7) and (8) are estimated by Monte Carlo integration. For a quantity f(r) and for the cell ξ the integral d^{3}r is approximated by (10)where we made explicit use of the mesh coordinates and where (ρ_{k},z_{k}) with k = 1,··· ,N are chosen randomly and uniformly within the cell domain.
Fig. 2 Quadtree mesh for a disc configuration. The disc parameters are those of model (b) described in Sect. 5.2 (see also Tables 2 and 4). The mesh refinement parameter η (Eqs. (7) and (8)) has been set to the high value 10^{3} in order to obtain a coarse mesh more easily represented. 

Open with DEXTER 
3.2. Symmetries
We can make use of the CSE symmetries to reduce the computation domain of an intensity map from Eq. (4) to only a fourth of it and consequently reduce the computation time.
Recalling the definition of a ray (Eq. (2)), we have two noticeable identities for any disc physical quantity Φ (e.g. , , T, n, ...) depending on rwhere Eq. (11) expresses the disc symmetry with respect to the y − z plane and Eq. (12) the point symmetry with respect to the origin of the model system of coordinates.
From the above identities it is straightforward to deduce their counterpart for the intensity map where . Note that the exponential factor in Eq. (14) has to be evaluated when computing I_{λ}(X,Y,i) anyway; no extra effort is required to derive I_{λ}(X,Y,i) from I_{λ}(X, − Y,i) except for the multiplication of I_{λ}(X, − Y,i) by this factor.
3.3. Intensity map
The fictitious image plane is split into a set of pixels whose positions X_{j} and Y_{k} are given by where Δ_{X} = Δ_{Y} is the pixel size in X and Y, and N is the number of pixels in X and Y, and where (17)where δ = −1 for N even and δ = 0 otherwise and “ ÷ ” stands for the integer division. Taking into account the symmetries mentioned in Sect. 3.2 only a fourth of the pixels need to be considered.
The evaluation of the integral in Eq. (4) is carried out for each pixel (X_{j},Y_{k}) and along the ray r_{s}(X_{j},Y_{k},i). The intersection points of the ray with the cell boundaries corresponds to a set of distances along the ray defined as where n_{cells} is the number of cells encountered along the ray, and Δs_{l} the distance crossed within the lth cell.
We estimate numerically the optical depth τ_{λ}(X,Y,i;s), defined in Eq. (3), via the midpoint rule quadrature by (20)where we defined for l = 0,··· ,n_{cells} − 1.
The numerical estimate of I_{λ}(X_{j},Y_{k},i) is obtained by (21)The results for all j and k can then be obtained from the discrete counterpart of the symmetry relations (13) and (14)
3.4. Interferometric observables
From the numerical estimate of I_{λ}(X_{j},Y_{k},i) given above we obtain (similarly to Eqs. (5) and (6)) the numerical fluxes and visibilities, which can be directly compared to the observed data. The numerical estimate of these quantities is again obtained through the midpoint rule.
The numerical flux F_{λ}(i) is computed by (24)The complex visibility is approximated numerically by (25)where B = (BcosΔ,BsinΔ) and R_{kl} = (X_{k},Y_{l}).
3.5. Artificial data generation
The procedure described below aims to mimic the observables of the VLTI/MIDI instrument: the flux F_{λ} (Eq. (24)) and the modulus of the visibility V_{λ} (Eq. (25)). The wavelengths and baselines chosen for the artificial data generation correspond to accessible values to VLTI/MIDI with the Unit Telescopes (UTs): λ_{j} = 7,8,9,10,11,12, and 13 μm (j = 1,··· ,n_{λ}; n_{λ} = 7), and (B_{k},PA_{k}) as shown in Table 1 (k = 1,··· ,n_{B}; n_{B} = 18). These values amount to 126 points covering the uvplane.
For a given intensity map at λ_{j}, F_{λj} and V_{λj} are taken as the expectation values of the simulated data. The observed flux is then generated assuming a Gaussian noise with an RMS (root mean square) corresponding to 10% relative error σ_{F}(j) = 0.1 × F_{λj}.
The artificial observed visibility amplitudes are obtained as (26)where ΔV_{k} is a wavelength independent shift that mimics the error in the observed visibilities, introduced by the calibration procedure commonly used in optical/IR interferometry. For each (B_{k},PA_{k}), ΔV_{k} is computed assuming a Gaussian noise with an RMS corresponding to 10% relative error (typical for VLTI/MIDI) σ_{V}(k) = 0.1 × ⟨V_{λ}(B_{k},PA_{k})⟩ , where ⟨V_{λ}(B_{k},PA_{k})⟩ is the wavelength mean visibility modulus.
Projected baselines. These values correspond to the baselines accessible from pairs of Unit Telescopes (UT) at ESOVLTI.
3.6. Model fitting and error estimate
We describe here the procedure adopted in order to simultaneously fit observed fluxes and visibilities using FRACS models defined by a given set of input parameters. This procedure is applied to artificial data in the next sections.
In order to quantify the discrepancy between the artificial observations ( and ) and the visibilities and fluxes from a given model (V_{λj}(B_{k},PA_{k}) and F_{λj}) we use the χ^{2} like quantities (27)and (28)To take into account both the midIR flux and the visibilities on the same level in the fitting process, we minimise the following sum (29)In the discussion below about the parameter and error determination we use the reduced χ^{2} defined by (for n_{free} free parameters).
From a minimising algorithm the bestfit model parameters can then be found by determining the minimum : . The “error” estimate is obtained from a thorough exploration of the parameter space volume, defined by a contour level , where has been chosen equal to 1. This volume can be interpreted as a confidence region. The quantity defined in Eq. (29) is a weighted sum of χ^{2} variables whose cumulative distribution function can be approximated by a gamma distribution (see Feiveson & Delaney 1968) with the same mean and variance. It is then possible to obtain a rough estimate of the confidence level associated with the confidence region given approximately by ≃2σ.
The size of the confidence region is determined by considering all possible pairs of parameters for a given fitted model and computing maps for each. The procedure to estimate the errors can be summarised as follows:

For a given map, i.e. for a given couple of parameters among the n_{free} × (n_{free} − 1)/2 possibilities, we identify the region bounded by around the minimum of this particular map.

The boundaries of the projection of these regions on each of the two parameter axis considered are recorded for each map.

The final errors on a given parameter are taken as the highest boundary values of the projected regions over all maps.
4. Astronomical test case: sgB[e] stars
In the following sections we apply FRACS to a theoretical interferometric study of dusty CSE of B[e] supergiants (sgB[e] in the nomenclature of Lamers et al. 1998). However, we emphasise that FRACS is in no way restricted to this particular class of objects.
sgB[e] stars reveal in particular a strong near or midIR excess caused by hot dust emission. There is evidence (e.g. Zickgraf et al. 1985) that the stellar environment, and in particular dust, could be confined within a circumstellar disc. Our purpose is to characterise this class of objects and derive not only geometrical parameters (e.g. inner dust radius, disc position angle and inclination) but also physical parameters such as temperature gradients, dust formation region, material density, ...
The physical description of the CSE chosen for our study is the wind model with equatorial density enhancement. This is a classical CSE model commonly adopted for sgB[e] (e.g. Porter 2003).
In order to compute the model intensity maps we need to parameterise the emissivity of the disc. Consistently with FRACS assumptions, we consider only dust thermal emission without scattering by dust grains and the gas contribution to the medium emissivity. In the rest of this section we characterise the emissivity by describing the dust density law, the absorption cross section, and the temperature structure of the CSE.
4.1. Mass loss and dust density
Dust is confined between the inner and outer radius R_{in} and R_{out} respectively. We assume a stationary and radial mass loss; physical quantities will consequently depend only on the radial component r and the colatitude θ. The disc symmetry axis coincides with the z axis of the model cartesian system of coordinates. The mass loss rate and velocity parametrisations are simplifications of the one adopted by Carciofi et al. (2010), and we refer the reader to their work for a complete description (see also Stee et al. 1995, for a similar description).
The mass loss rate per unit solid angle, at colatitude θ, is parameterised as follows (30)with the help of two dimensionless parameters A_{1} and m.
Even though our computations make no explicit use of the radial velocity field v_{r}(θ) (assumed to have reached the terminal velocity v_{∞}(θ) in the region under considerations, i.e. v_{r}(θ) ≈ v_{∞}(θ)), the dust density depends on v_{r}(θ) parameterised in a similar fashion (31)where we have introduced the supplementary dimensionless parameters A_{2}. From Eqs. (30) and (31) we see that A_{1} and A_{2} are the relative differences of the values of and v_{r}(θ) at the equator and the pole (relatively to the pole).
From the mass continuity equation one obtains the number density of dust grains (32)where n_{in} is the dust grain number density at R_{in} in the disc equatorial plane. In Eq. (32), the parameter m controls how fast the density drops from the equator to the pole, defining an equatorial density enhancement (disclike structure).
Consistent with the accepted conditions for dust formation (Carciofi et al. 2010; Porter 2003) we assume that the dust can survive only in the denser parts of the disc. We thus define a dusty disc opening angle Δθ_{d} determined by the latitudes for which the mass loss rate has dropped to half of its equatorial value: (33)To summarise, the dust grains only exist (i.e., n(r,θ) ≠ 0) in the regions bounded by R_{in} ≤ r ≤ R_{out} and by .
4.2. Dust opacities
The absorption cross section for the dust grains is obtained from the Mie (1908) theory. The Mie absorption cross sections are computed from the optical indices of astronomical silicate (Draine & Lee 1984). Note that since scattering is neglected, , with being the extinction cross section.
For a powerlaw size distribution function according to Mathis et al. (1977) the mean cross sections (e.g. for ) are given by (34)where a_{min} and a_{max} are the minimum and maximum radii for the dust grains under consideration and β is the exponent of the powerlaw. The computation of the cross section in Eq. (34) was performed with the help of the Wiscombe (1980) algorithm.
4.3. Temperature structure
The dust temperature is assumed to be unique (i.e. independent of grain size) and described by a powerlaw (35)where T_{in} is the temperature at the disc inner radius R_{in}. We note that γ is not necessarily a free parameter because in the optically thin regime (large wavelength and radius) the temperature goes as with δ ≃ 1 (see Lamers & Cassinelli 1999).
4.4. Central continuum emission
The continuum emission from the central regions is composed by the emission from the star and from the close ionised gas (freefree and freebound emission). This central source emission is confined to a small region of radius R_{s} (≪ R_{in}), which is unresolved (angular sizes of a few milliarcseconds) by midIR interferometers. Thus, in our modelling R_{s} is simply a scaling factor of the problem fixed to a typical radius value for massive stars. The specific intensity (in W m^{2} μm^{1} str^{1}) of this central source is parameterised as follows (36)where is the specific intensity at a reference wavelength λ_{0} (= 10 μm in the following), and α gives the spectral dependence of the continuum radiation. In the midIR its value is expected to lie between α = 4 (pure black body) and α ≃ 2.6 (freefree emission) for an electron density proportional to r^{2} (Panagia & Felli 1975; Felli & Panagia 1981).
5. Study of the tested models
Model parameters. This table lists the parameters of 12 different models.
Following the description in the last sections we describe here the chosen sgB[e] model parameters used to simulate VLTI/MIDI observations (visibilities and fluxes) and the corresponding analysis, i.e. model fitting, using FRACS. The list of chosen parameters is summarised in Table 2. Two types of numerical tests are presented. Firstly, synthetic midIR interferometric data are generated from FRACS itself. In that way, it is possible to determine what information the midIR interferometric data contain under the optimistic assumption that we do have the true model. Secondly, this study is supplemented by the comparison of FRACS to a Monte Carlo radiative transfer computation. This confirms that FRACS can indeed mimic, under appropriate conditions, the results of a more sophisticated code as seen from the midIR interferometric eye.
5.1. Parameter description
The distance to the simulated object has been fixed to d = 1 pc, which is a typical distance for Galactic sgB[e].
The inner radius R_{in} = 30 R_{s} = 1800 R_{⊙} value was chosen by considering the location of the hottest dust grains (see Lamers & Cassinelli 1999) with a condensation temperature of 1500 K assumed to be the T_{in} value. The value of R_{out} cannot be determined from the midIR data and has been fixed to 3000 R_{s} = 1.8 × 10^{5} R_{⊙}. The temperature gradient γ was fixed to 0.75 according to Porter (2003). PA_{d} was fixed to 125°.
The central source emission is supposed to have a radius R_{s} = 60 R_{⊙}. We recall that the central region is unresolved by the interferometer and that its radiation describes both the stellar and inner gas contribution to the continuum midIR emission. The specific intensity of this central source has been chosen to be 6500 W m^{2} μm^{1} str^{1}. If the central source was a pure blackbody this value would correspond to the 10 μm emission of a blackbody with an effective temperature around ≃8000 K. However, this central emission is not a pure blackbody, and we adopt the spectral dependence of the central source emission to be α = 3, which is a compromise between α = 4 for a pure blackbody and a value of ≃2.6 for freefree emission (Panagia & Felli 1975; Felli & Panagia 1981).
Spectroscopic observations of Hα and forbidden line emissions from B[e] CSE (Zickgraf 2003) reveal that typical values for A_{2} are expected to range from −0.95 to −0.75. We adopt the value −0.8 in our models. According to Lamers & Waters (1987), the values of A_{1} range from 10^{2} to 10^{4} in most cases (though values as low as 10 are not excluded). With this high value of A_{1} the factor 1 + A_{1}(sinθ)^{1/m} in Eq. (32) of n(r,θ) is approximatively given by A_{1}(sinθ)^{1/m} for all pertinent values of θ, i.e. those close to π/2 within the disc. This leads to an evident degeneracy in n_{in} × A_{1} in n(r,θ): we are only sensible to the product of the two parameters as a scaling factor for the density. Therefore, the value of A_{1} is assumed to be fixed to 150.
Relative errors (given in %) on the parameters for models (a) et (b) (see text for description of models).
Constraints on the model parameters.
To define the dust opacities the chosen value for β is that of Mathis et al. (1977), i.e. β = −3.5. Because some sgB[e] show weak 9.7 μm silicate features in their spectrum (e.g. Porter 2003; Domiciano de Souza et al. 2007) we chose to use large grains in our test models: a_{min} = 0.5 μm and a_{max} = 50 μm. However, with this particular choice of large grains, the average albedo from 7 to 13 μm is 6.4%, with the highest value reached at 7 μm. We have checked with a Monte Carlo (MC) simulations (see Sect. 5.3) that the effect of scattering on our primary observables, visibilities, and fluxes is indeed negligeable by comparing the results obtained by switching the scattering process off and on^{2}. The mean relative differences are 3.5% and 3.0% for the visibilities and the fluxes respectively. These values must be compared to the effect of random noise in the MC simulation, estimated to be of the same order and to experimental errors, typically ~10% for the visibilities and fluxes. We underline that whenever the albedo can be neglected, it is theoretically safe to compute visibilities and fluxes from the consideration presented in Sect. 2, in any other situations the effect of scattering on the observable must be carefully tested.
The parameters n_{in}, m and i were set to different values defining 12 test models to be analysed from their corresponding simulated data. Two n_{in} values (0.015 m^{3} and 0.15 m^{3}) have been chosen in order to have an approximate discdust optical depth in the equatorial plane (from R_{in} to R_{out}) close to ≃0.1 and ≃1 in the wavelength range considered (from 7 μm to 13 μm). These values corresponds to a mass loss rate of Ṁ ≃ 2.5 × 10^{ −7··· −6} M_{⊙} yr^{1}. Two m values were chosen corresponding to a wide and a narrow opening angle, i.e. Δθ_{d} = 10° and Δθ_{d} = 60°. Three inclinations i were tested (20°, 50°, and 90°) corresponding to discs seen close to poleon, intermediate inclination, and equatoron. These values of n_{in}, m, i, together with the parameters fixed above, define 12 test models that will be studied below.
From these 12 test models we have generated 12 sets of artificial VLTI/MIDI observations (visibilities and fluxes) following the procedure described in Sect. 3.5. We do not aim to present an exhaustive revue of all types of sgB[e] CSE. Rather, we focussed on the analysis of the parameter constraints one can hope to obtain from present and nearfuture midIR spectrointerferometry. The quantitative estimate of these constraints is derived from a systematic analysis of the variations with the parameters.
In our model fitting and analysis we concentrate on 10 free parameters (n_{free} = 10) that can be set into four different groups:

the geometrical parameters: PA_{d}, i and R_{in};

the parameters related to the central source: and α;

those describing the temperature structure: T_{in} and γ;

and the number density of dust grains: A_{2}, Δθ_{d} (or equivalently m) and n_{in}.
The remaining parameters of the model are in general loosely constrained by midIR interferometric observations so that we kept them fixed to the values described above.
Fig. 3 Simulated VLTI/MIDI midIR object flux. The numerically generated data (left model a, right model b) are shown as circles with error bars. The values of the fluxes for the bestfit models are represented as crosses. 

Open with DEXTER 
Fig. 4 Visibilities of the artificial sgB[e] circumstellar environment (model a). The visibility variations with the wavelength are shown for each baseline specified by the value of the projected baseline and the position angle on the sky. The circles represent the simulated observations, and the solid curves represent the bestfit model. 

Open with DEXTER 
Fig. 5 Visibilities of the artificial sgB[e] circumstellar environment (model b). 

Open with DEXTER 
5.2. Model fitting and analysis of the 12 test models
We describe here the data analysis procedure adopted to study our 12 test models. The results of our analysis are summarised in Tables 3 and 4, and their physical interpretation is presented in Sect. 6.
As a first step we chose 2 of the 12 models, hereafter called models (a) and (b), to be exhaustively studied from a complete model fitting procedure. As an example we show the simulated observed midIR fluxes and visibilities for model (a) and (b) in Figs. 3−5. The parameters of models (a) and (b) are those of Table 2 with Δθ_{d} = 60°, i = 50° and n_{in} fixed to the values 0.015 m^{3} and 0.15 m^{3} respectively. These two models are those presenting some of the best constrained model parameters for the dust CSE. On the other hand, the contribution of the central regions to the total flux and visibilities is quite different in models (a) and (b) (see discussion in Sect. 6).
The study of models (a) and (b) have thus been performed as for real interferometric observations. The bestfit values of the parameters have been obtained by the LevenbergMarquardt algorithm with a stopping criterion corresponding to a relative decrease in of 10^{3}.
The errors on each model parameter have been obtained following the methods described in Sect. 3.6. The maps have been computed with a resolution of 21 × 21 around the bestfit values of the parameters. The map sizes have been adjusted in order to enclose the contour. This adjustment was performed until an upper limit for the map size of 100% of the bestfit parameter values was reached. This amounts to the computation of 3.969 × 10^{4} different models. The results, namely the mean relative error up to 100%, for these two particular models are summarised in Table 3.
The other ten models (numbered from 1 to 10 in Table 4) have been used in order to get some quantitative (but limited) information about how the uncertainties of the fitted parameters evolve as a function of three disc characteristics: its optical depth (τ by means of n_{in} parameter), its inclination (i) and its opening angle (Δθ_{d}, controlled by m). To perform this study we have decided to limit the exploration of the space parameter in a relative range of 25% on both sides from the model parameters. In order to reduce the computation time, the maps were not generated around the bestfit parameters which would have required to compute several thousands models more but around the true parameters themselves. This procedure has the supplementary advantage that we do not rely on any specific minimisation algorithm. We checked that estimating the bestfit parameters from the true ones is reasonable within a few percents using the LevenbergMarquardt algorithm with a stopping criterion corresponding to a relative decrease in of 10^{3}. The resolution of the maps have been reduced to 15 × 15. The total number of models to be computed is as large as 1.0125 × 10^{5}.
5.3. Comparison with a Monte Carlo simulation
bestfitting FRACS parameters from artificial data generated with a Monte Carlo code.
We generated synthetic data with the help of a Monte Carlo (MC) radiative transfer code (Niccolini & Alcolea 2006) for model b (see above) for the seven wavelengths considered in the problem and the baselines of Table 1. Again, the adopted procedure to generate the midIR interferometric data follows the considerations of Sect. 3.5. In the MC code, the source of photons is described by a blackbody sphere of radius R_{s} = 60 R_{⊙} and an effective temperature of T_{eff} = 8000 K. The temperature of the CSE is not prescribed but computed from the Lucy (1999) mean intensity estimator. This choice of T_{eff} gives at the inner radius of model b a dust temperature of ≃1150 K lower than the sublimation temperature. In this way, we can test if in the fitting process using FRACS, a spurious effect might not lead the minimisation algorithm to reach the upper limit for T_{in} of 1500 K, corresponding to the adopted dust sublimation temperature.
We obtained the best fitting parameters for the CSE model described in Sect. 4 with FRACS. For a comparison with the MC code, α has been set and fixed to 4 corresponding to the value of a blackbody. Depending of the disc optical depth, the temperature structure may show two separate regimes corresponding to (1) the inner regions with the strongest temperature gradient, optically thick to the stellar radiation and (2) the outer regions optically thin to the disc radiation with a flatter temperature gradient. In order to determine if midIR interferometric data are sensitive to two temperature regimes, we tested the effect of two parameterisations of the temperature structure: the unique powerlaw of Eq. (35) and a generalisation to two powerlaws with a transition radius, R_{T}, and a second exponent γ′(37)for r ≥ R_{T}.
The bestfitting parameters for both parameterisations are shown in Table 5. The images of the disc at 10 μm generated with the MC code and their corresponding FRACS counterpart (bestfitting model) are shown in Fig. 6 for comparison.
Fig. 6 Disc images at 10 μm. a) Image computed with the help of the Monte Carlo radiative transfer code. b) Image of the bestfitting model with two powerlaws (parameters of the fourth column in Table 5) obtained with FRACS. 

Open with DEXTER 
6. Discussion
We first discuss the uncertainties in the parameters derived for the 12 models studied in Sect. 5.2. For each model we divided the parameters into three groups associated to a given level of constraints expressed by the relative errors: below 10%, between 10% and 25%, and above 25%. This information is summarised in Table 4. The exact relative errors for the two models studied in detail (models a and b) are shown in Table 3. Then, we analyse the results of Sect. 5.3 obtained from the bestfit of the data simulated with the MC radiative transfer code.
6.1. Central source
Table 4 shows that the central source parameters ( and α) can only be constrained with relative uncertainties ≥10% for all test models. A deeper and more quantitative investigation of these parameters can be obtained from models (a) and (b). From Table 3 it can be seen that and α are much better constrained for model (a) (27% and 13%, resp.) than for model (b)(relative errors ≃100%).
The key quantity for a good constraint for the central source parameters ( and α) is simply the relative contribution of the flux of the central source to the total flux of the object (source and disc). Indeed, the models in Table 3 only differ by this relative flux contribution of 5.3% in model (a) (τ = 0.1), while it is only 0.7% in model (b) (τ = 1).
Our analysis thus shows that interferometric data can constrain and α with a relative precision of ≃15%−30% even when the central source contributes to (only) a few percent of the total midIR flux.
6.2. Geometrical parameters
The parameters PA_{d}, i, and R_{in} are those usually estimated from simple geometrical models (e.g. ellipses, Gaussians). However, their determination from geometrical models is quite limited, in particular for i, for which only an estimate can be derived from the axisratio of an ellipse, for example. In addition, the estimate of i from a simple analytical model such as a flat ellipse is only valid for configurations far from the equator (intermediate to low i). The use of a more physical and geometrically consistent model such as FRACS allows us to relax this constraint and makes the determination of i possible for all viewing configurations.
As expected, PA_{d} and i are better determined if the inclination of the disc with respect to the line of sight is away from poleon (high i). In Fig. 8 we can clearly see this behaviour from the maps involving PA_{d} and i. Moreover, the uncertainties on PA_{d} and i do not seem to be strongly dependent on τ (equivalently n_{in}) and Δθ_{d} (equivalently m) for all models.
The inner dust radius R_{in} is not strongly dependent on any parameter (τ, Δθ_{d} or i), being very well constrained (better than 10%) for most tested models.
6.3. Temperature
The parameters related to the temperature structure of the CSE, T_{in} and γ, are well constrained in most models, with relative errors below 20% and 12% for both models (a) and (b). Indeed, γ has a strong impact on the IR emission across the disc, and consequently this parameter has a direct influence on the visibilities (see Figs. 4, 5). T_{in} has a lower influence, compared to γ, on the shape of the monochromatic image (radial dependence of intensity) and can be mainly considered as a scaling factor to it. On the other hand, the midIR flux imposes stronger constraints on T_{in}. From Table 4 we see that the CSE’s temperature structure is not highly dependent on τ (n_{in}) and Δθ_{d} for tested models.
6.4. Number density of dust grains
The parameters related to the density law, that is to say m, n_{in} and A_{2}, seems to be rather poorly constrained from the midIR data alone. From the results of model (a) and (b) corresponding to an intermediary inclination i = 50°, we found that only n_{in} is constrained somewhat moderately with a mean relative error of 46%. For m and A_{2}, according to the results of Table 3 it seems that nevertheless, upper limits to their values can be determined. Note that because A_{2} is bounded (−1 ≤ A_{2} ≤ 0) the mean relative errors, 53% and 59% for model (a) and (b) respectively, correspond approximatively to the limit values of A_{2}, which is consequently not constrained.
Table 4 confirms this trend for m, n_{in} and A_{2} at least for the situations explored via the models presented here. From all maps computed within ±25% of the true value, we always found that the mean relative error to these parameters is larger than 25% with no hint that it could be close to these limits.
From Fig. 9, comparing the maps for all pairs of n_{in}, m and A_{2}, we can see that the contours get sharper around the minimum value for model (a) (corresponding to lower optical depths along the line of sight) than for model (b). Indeed, the constraints on n_{in} and m are improved for lower optical depths, or equivalently for lower disc masses. Indeed, when the disc mass (or optical depth) decreases, the flux (midIR flux, intensity maps) emitted by the disc reflects the mass of the disc, while for high optical depths we only probe the regions of the disc very close to the projected surface revealed to the observer. A_{2}, however, is unaffected by the change in disc mass and remains undetermined anyway. From Fig. 9 it can be seen that n_{in}, m and A_{2} are strongly correlated. This is expected from the expression of the density (see Eq. (32)) depending on these parameters. However, this dependence and the final correlation between these parameters are related through the computation of the visibilities and the midIR flux, as well as the comparison to the data and is, therefore, not straightforward.
To improve the situation concerning n_{in}, m and A_{2}, the midIR data can be supplemented by other types of observations such as for instance spectroscopic data, from which one can better determine A_{2} (e.g. see Chesneau et al. 2005). We tested the effect of fixing the value of A_{2}, or equivalently of assuming that A_{2} is fully determined, in the process of estimating the errors of the other parameters. For model (a), the relative errors on n_{in} and m go down to 33% and 71% respectively while for model (b), n_{in} and m are determined with an accuracy reaching 95% and 78% respectively. The precision to which other parameters are determined is not affected by the determination of A_{2}.
We also tested the influence of the determination of n_{in}, m and A_{2} on other parameters by fixing their values and estimatng the relative errors on the remaining parameters for model (a) and (b). Only and T_{in} are more strongly affected by the determination of n_{in}, m and A_{2}: for model (a) (resp. model b) gets determined down to 19% (resp. 80%) and T_{in} down to 9% (resp. 18%). The influence is stronger with lower disc mass (model a compared to model b). This effect can be explained because if we have a good determination of the disc mass because we know n_{in}, m, and A_{2}, the determination of the parameters that scale the source and disc fluxes is improved accordingly for the visibility and the midIR flux.
n_{in}, m and A_{2} shape the density structure of the circumstellar medium. Though they are not well constrained, they certainly have a strong influence on the temperature structure, which in turn is very well constrained. For the particular case of sgB[e] circumstellar discs, a natural evolution of FRACS is to include the direct heating of the medium by the central source of radiation assuming that the disc is optically thin to its own radiation. The temperature structure would not be parameterised, and its good determination would certainly put better constraints on n_{in}, m and A_{2}, while keeping an affordable computation time for the modelfitting procedure. This will be the purpose of a subsequent work.
Finally, one can derive a ranking of the parameter constraints according to two criteria: first the parameter must be constrained within the prescribed limits (100% for model a and b and 25% for model 1 to 10) and second the mean relative error must be as low as possible. The bestfitted parameters, most of the time according to these criteraria are by decreasing order of best determination: R_{in}, PA_{d}, γ, T_{in}, i, α, , n_{in}, m and A_{2}. This tendency can be seen in Table 4.
6.5. bestfit to the MC simulation
Fig. 7 Temperature of the CSE. The solid line represents the best fit (last column in Table 5) with a unique powerlaw, the dashed line the bestfit with two powerlaws (fourth column in Table 5) and the dotdashed line the MC results. The shaded region represents the possible domain for a unique powerlaw by taking into account the errors estimated in Table 3. 

Open with DEXTER 
The values obtained for the two types of temperature parametrisations (one and two powerlaws) are quite similar: 0.79 and 0.73 respectively. Regarding the data, both temperature parametrisations are indeed acceptable. In addition, these results show that we can actually obtain very good fits from data sets based on more physically consistent scenarios. A complete error analysis and study of the parameter determination has been presented in the previous sections for data generated from FRACS and will not be repeated here. In particular, parameter confidence intervals, from which errors were derived, have already been estimated. Here, we will instead focus on the true errors, i.e. the differences between the true model parameters and the bestfitting values for the parameters (see Table 5). The two types of errors must not be confounded. The true errors reflect the capability of FRACS to mimic the midIR interferometric data regarding the information it provides. Of course, with the sparse uvplane coverage inherent to this kind of data as well as the experimental noise, one should not expect a full agreement of the fitted and the true parameters: they are indeed different.
From Table 5, we see that the geometrical parameters, PA_{d}, i and R_{in}, can be almost exactly recovered as expected. The source specific intensity , and the parameters related to the density, A_{2}, n_{in} and m can be recovered fairly well and have bestfitting values close to the true parameters.
Fig. 8 Evolution of m and PA_{d} with the inclination i. Left: maps for the couple m and PA_{d}; right: χ^{2} maps for the couple m and i. Contours are drawn for , with . From top to bottom the inclination i takes the value 20°, 50° and 90°. The results correspond to model 4, 6, and 10. The limits of the maps have been set to ±25% of the true values of the parameters. 

Open with DEXTER 
Fig. 9 maps for the parameters n_{in}, A_{2} and m. The results presented here are those of model (a) (τ = 0.1, left part) and model (b) (τ = 1, right part). Contours are drawn for , with . The three possible maps corresponding to the combination of these parameters are represented. These three parameters are better constrained in model (a). 

Open with DEXTER 
Fitting the temperature computed with the MC code with a simple powerlaw, we obtain γ ≃ 0.64. This value is close to those of the best fitting models, especially with a unique powerlaw (6% relative difference). For comparison, the actual mean temperature gradient as derived from the MC simulation is ≃0.60. For this particular data set, the values of γ′ and R_{T} recovered by FRACS differ by 28% and 45% respectively from the actual values. This again confirms the sensibility of the interferometric data to the temperature structure mostly in the inner (r ≲ R_{T}) regions of the disc. The bestfitting models (fourth and last columns in Table 5) as well as the MC results are shown in Fig. 7. Regarding the errors (estimated from the results given in Table 3) shown as a shaded area, we can see that both temperature parameterisations are essentially the same and show a better agreement with the MC results in the inner than in the outer regions of the disc.
We considered a “truncated” model with two powerlaws (with parameter values listed in the third column of Table 5) in which the CSEs emission for r ≥ R_{T}, the “outer” regions, has been set to 0. We then compared the visibilities and the fluxes of this truncated model to the same model including the outer region emission. We obtained relative differences, averaged over all considered wavelengths and baselines (Table 1), of 18% and 17% respectively. These relative differences are larger, but are still close to the noise level. For this reason, one cannot expect to obtain much information on the outer temperature gradient γ′, at least for the particular configuration we considered.
7. Conclusion
We proposed and described here a new numerical tool to interpret midIR interferometric data. Even though we focussed on the special case of circumstellar disc observations, the numerical techniques have been developed with the aim to be as general as possible. The methods we employ rely on both parameterised physical models and the raytracing technique. The need for such a tool is evident because the nature of interferometric data imposes an interpretation through a model of the object to obtain any kind of information. On one hand, MonteCarlo radiative transfer methods require too much computation time to associate the modelfitting to an automatic minimum search method. On the other hand, purely geometrical function fitting (such as ellipses or Gaussians) are too simple to envisage to obtain physical constraints on the observed disc. Hence, a tool like FRACS fills a blank in the model fitting approach for midIR interferometric data interpretation. The main advantages of FRACS are its speed and its flexibility, allowing us to test different physical models. Moreover, an exploration of the parameter space can be performed in different manners and can lead to an estimate of the sensitivity of the fit to the different model parameters, i.e. a realistic error estimate.
We applied these techniques to the special astrophysical case of B[e] star circumstellar environments by generating artificial data in order to analyse beforehand what constraints can be obtained on each parameter of the particular disc model in this work. The techniques will then be applied to real interferometric data of a sgB[e] CSE in a sequel to this paper.
We showed in our analysis that the “geometrical” parameters such as R_{in}, PA_{d} and i can be determined with an accuracy ≲15%. MidIR interferometric data give access to a mean temperature gradient: the temperature structure (T_{in} and γ) can be very well determined (within ≲20% and ≲10% respectively). It is possible to have access to the central source emission (with an accuracy ≳30%) when it has a significant contribution to the total flux of the object (a few % are sufficient). The remaining parameters of our disc model, namely n_{in}, m and A_{2} are not very well constrained by MIDI data alone. n_{in} is at best determined with an accuracy of about ≳50% in some cases. If A_{2} can be estimated through spectroscopic observations, then the picture about the n_{in} and m determination improves somewhat.
FRACS can be used mainly for two purposes. First, it can be used by itself to try and determine physical quantities of the circumstellar matter. Admittedly, it is not a selfconsistent model, i.e. the radiative transfer is not solved because the temperature structure is parameterised. From the usual habits in the interpretation of interferometric data it is nevertheless a step beyond the commonly use of toy models or very simple analytical models.
This approach has indeed been very successful in the millimetric wavelength range (e.g. see Guilloteau & Dutrey 1998). Second, it can be viewed as a mean to prepare the work of data fitting with a more elaborate model (such as a Monte Carlo radiative transfer code for instance) and to provide a good starting point.
FRACS is a tool that can help in the process of interpreting and/or preparing observations with secondgeneration VLTI instruments such as the MultiAperTure midInfrared SpectroScopic Experiment (MATISSE) project (Lopez et al. 2006). In this respect, FRACS is not restricted to the midIR, and submillimeter interferometric data obtained with the Atacama Large Millimeter Array (ALMA) for instance can be tackled.
Acknowledgments
We thank the anonymous referee for her/his constructive comments that lead to significant improvements of the manuscript. We also would like to thank Alex Carciofi, Olga Suárez Fernández, Andrei and Ivan Belokogne for fruitful discussions and careful proof reading.
The Monte Carlo simulations have been carried out on a computer financed by the BQR grant of the Observatoire de la Côte d’Azur.
This work is dedicated to Lucien.
References
 Bianchi, S. 2008, A&A, 490, 461 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Carciofi, A. C., Miroshnichenko, A. S., & Bjorkman, J. E. 2010, ApJ, 721, 1079 [NASA ADS] [CrossRef] [Google Scholar]
 Chesneau, O., Meilland, A., Rivinius, T., et al. 2005, A&A, 435, 275 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 di Folco, E., Dutrey, A., Chesneau, O., et al. 2009, A&A, 500, 1065 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Domiciano de Souza, A., Driebe, T., Chesneau, O., et al. 2007, A&A, 464, 81 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Draine, B. T., & Lee, H. M. 1984, ApJ, 285, 89 [NASA ADS] [CrossRef] [Google Scholar]
 Feiveson, A. H., & Delaney, F. C. 1968, The distribution and properties of a weighted sum of chi squares, Tech. rep., National Aeronautics and Space Administration [Google Scholar]
 Felli, M., & Panagia, N. 1981, A&A, 102, 424 [NASA ADS] [Google Scholar]
 Guilloteau, S., & Dutrey, A. 1998, A&A, 339, 467 [NASA ADS] [Google Scholar]
 Jonsson, P. 2006, MNRAS, 372, 2 [NASA ADS] [CrossRef] [Google Scholar]
 Kurosawa, R., & Hillier, D. J. 2001, A&A, 379, 336 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lachaume, R., Preibisch, T., Driebe, T., & Weigelt, G. 2007, A&A, 469, 587 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Lamers, H. J. G. L. M., & Cassinelli, J. P. 1999, Introduction to Stellar Winds, ed. J. P. Lamers, & H. J. G. L. M. Cassinelli [Google Scholar]
 Lamers, H. J. G. L. M., & Waters, L. B. F. M. 1987, A&A, 182, 80 [NASA ADS] [Google Scholar]
 Lamers, H. J. G. L. M., Zickgraf, F., de Winter, D., Houziaux, L., & Zorec, J. 1998, A&A, 340, 117 [NASA ADS] [Google Scholar]
 Leinert, C., Graser, U., Przygodda, F., et al. 2003, Ap&SS, 286, 73 [NASA ADS] [CrossRef] [Google Scholar]
 Leinert, C., van Boekel, R., Waters, L. B. F. M., et al. 2004, A&A, 423, 537 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Levenberg, K. 1944, The Quarterly of Applied Mathematics, 2, 164 [Google Scholar]
 Lopez, B., Wolf, S., Lagarde, S., et al. 2006, in SPIE Conf. Ser., 6268 [Google Scholar]
 Lucy, L. B. 1999, A&A, 344, 282 [NASA ADS] [Google Scholar]
 Malbet, F., Lachaume, R., Berger, J., et al. 2005, A&A, 437, 627 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Marquardt, D. 1963, SIAM J. Appl. Math., 11, 431 [CrossRef] [MathSciNet] [Google Scholar]
 Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 217, 425 [NASA ADS] [CrossRef] [Google Scholar]
 Mie, G. 1908, Ann. Phys., 25, 377 [NASA ADS] [CrossRef] [Google Scholar]
 Niccolini, G., & Alcolea, J. 2006, A&A, 456, 1 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Ohnaka, K., Driebe, T., Hofmann, K., et al. 2006, A&A, 445, 1015 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Panagia, N., & Felli, M. 1975, A&A, 39, 1 [NASA ADS] [Google Scholar]
 Porter, J. M. 2003, A&A, 398, 631 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Press, W. H., Teukolsky, S. A., Vetterling, W. T., & Flannery, B. P. 1992, Numerical recipes in C. The art of scientific computing (Cambridge: University Press) [Google Scholar]
 Stee, P., de Araujo, F. X., Vakili, F., et al. 1995, A&A, 300, 219 [NASA ADS] [Google Scholar]
 Steinacker, J., Bacmann, A., & Henning, T. 2006, ApJ, 645, 920 [NASA ADS] [CrossRef] [Google Scholar]
 Wiscombe, W. J. 1980, Appl. Opt., 19, 1505 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
 Wolf, S., Henning, T., & Stecklum, B. 1999, A&A, 349, 839 [NASA ADS] [Google Scholar]
 Zickgraf, F. 2003, A&A, 408, 257 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Zickgraf, F., Wolf, B., Stahl, O., Leitherer, C., & Klare, G. 1985, A&A, 143, 421 [NASA ADS] [Google Scholar]
All Tables
Projected baselines. These values correspond to the baselines accessible from pairs of Unit Telescopes (UT) at ESOVLTI.
Relative errors (given in %) on the parameters for models (a) et (b) (see text for description of models).
bestfitting FRACS parameters from artificial data generated with a Monte Carlo code.
All Figures
Fig. 1 Coordinate systems. The shaded ellipse represents a disc viewed by the observer. 

Open with DEXTER  
In the text 
Fig. 2 Quadtree mesh for a disc configuration. The disc parameters are those of model (b) described in Sect. 5.2 (see also Tables 2 and 4). The mesh refinement parameter η (Eqs. (7) and (8)) has been set to the high value 10^{3} in order to obtain a coarse mesh more easily represented. 

Open with DEXTER  
In the text 
Fig. 3 Simulated VLTI/MIDI midIR object flux. The numerically generated data (left model a, right model b) are shown as circles with error bars. The values of the fluxes for the bestfit models are represented as crosses. 

Open with DEXTER  
In the text 
Fig. 4 Visibilities of the artificial sgB[e] circumstellar environment (model a). The visibility variations with the wavelength are shown for each baseline specified by the value of the projected baseline and the position angle on the sky. The circles represent the simulated observations, and the solid curves represent the bestfit model. 

Open with DEXTER  
In the text 
Fig. 5 Visibilities of the artificial sgB[e] circumstellar environment (model b). 

Open with DEXTER  
In the text 
Fig. 6 Disc images at 10 μm. a) Image computed with the help of the Monte Carlo radiative transfer code. b) Image of the bestfitting model with two powerlaws (parameters of the fourth column in Table 5) obtained with FRACS. 

Open with DEXTER  
In the text 
Fig. 7 Temperature of the CSE. The solid line represents the best fit (last column in Table 5) with a unique powerlaw, the dashed line the bestfit with two powerlaws (fourth column in Table 5) and the dotdashed line the MC results. The shaded region represents the possible domain for a unique powerlaw by taking into account the errors estimated in Table 3. 

Open with DEXTER  
In the text 
Fig. 8 Evolution of m and PA_{d} with the inclination i. Left: maps for the couple m and PA_{d}; right: χ^{2} maps for the couple m and i. Contours are drawn for , with . From top to bottom the inclination i takes the value 20°, 50° and 90°. The results correspond to model 4, 6, and 10. The limits of the maps have been set to ±25% of the true values of the parameters. 

Open with DEXTER  
In the text 
Fig. 9 maps for the parameters n_{in}, A_{2} and m. The results presented here are those of model (a) (τ = 0.1, left part) and model (b) (τ = 1, right part). Contours are drawn for , with . The three possible maps corresponding to the combination of these parameters are represented. These three parameters are better constrained in model (a). 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.