© P. Mayer et al. 2022

1. Introduction

The massive multiple system HD 93206 is a very unusual and rare object, which could potentially be very important for current theories of the origin of multiple systems and their evolution. If its basic physical properties can be obtained with sufficient accuracy, we may be able to estimate the evolutionary age of individual stars and to decide whether they have a common origin. It is also important to verify whether or not one of the stars has already had time to evolve to the stage of evolutionary expansion and Roche-lobe overflow in the binary subsystem. Our study is a new and rather detailed contribution to the efforts of a number of previous investigators to achieve these goals. We collected all available previous individual observations and obtained numerous new ones. Their analyses are presented here.

2. Current knowledge about HD 93206

The massive multiple system known as HD 93206 (HIP 52526), a member of the open cluster Collinder 228, consists of at least nine stars. We adopt for them the same notation as that used in the Washington Double Stars catalogue (WDS)¹. Their positions in the sky are shown in Fig. 1 and their magnitudes, if known, are listed in Table 1. It is clear that the flux from the object is dominated by the components of the quadruple system Aa1+Aa2 and Ac1+Ac2, all other known components being too faint to affect the photometric measurements, even if observed with an electronic detector with larger pixel size. The rest of this study will therefore deal with the close quadruple system only. A guideline for orientation with the different names of these four components as used by different investigators is provided in Table 2.

Fig. 1. Sketch of the multiple system HD 93206 drawn in arcsec scale on both axes. We note the expanded scale on the ordinate. The red circle is centred on the position of the unresolved Aa1-Aa2 spectroscopic binary. North is up, east is to the left.

Throughout this paper, we also use the abbreviated form for heliocentric Julian dates, RJD = HJD-2400000.0, to avoid confusion, which sometimes happens when the half-day shift introduced by the modified Julian dates is overlooked.

Periodic light variations of HD 93206 with a 60 period were discovered by Walker & Marino (1972) and interpreted as binary eclipses with a β Lyr-type light curve. The eclipsing binary is now known as QZ Car. Conti et al. (1977) noted the presence of double lines in one photographic spectrum. Subsequent spectroscopic studies carried out by Morrison & Conti (1979), Leung et al. (1979), and Morrison & Conti (1980) demonstrated that there are two sets of spectral lines, one showing the radial-velocity (RV) variations with the eclipsing-binary 60 period, and the other with a 207 period, indicative of the presence of another spectroscopic binary. Leung et al. (1979) also derived the first light-curve solution and concluded that the eclipsing binary is a semi-detached system. This view was adopted by Morrison & Conti (1980), who also obtained a series of yellow-red spectra covering the Hα line. These authors demonstrated the complexity of the observed spectra and problems of their analyses. Large width of spectral lines in combination with line blends causes differences in the observed systemic velocities. In addition, the velocity amplitudes of lower- and higher-ionisation lines were different. Morrison & Conti (1980) also found the presence of a sharp nebular emission from the Carina nebula at a RV from −27 to −28 km s⁻¹ in Hα and [O II] 3729 Å. They argued that at least part of the complicated Hα emission moves with the eclipsing Ac system.

Mayer et al. (2001) obtained new spectral and photometric observations and improved ephemerides of both orbits. Considering the light-time effect, these authors tentatively suggested that the orbital period of the Aa–Ac system can be about 50 yr. Stickland & Lloyd (2001) obtained RVs of components Aa1, Ac1, and for the first time also Ac2 from nine Short Wavelength Prime (SWP) International Ultraviolet Explorer (IUE) spectra and derived new orbital elements of both orbits. For the Ac1–Ac2 subsystem, they obtained a mass ratio M_Ac1/M_Ac2 = 1.07, and large values M_Ac1sin³i = 43.6 ± 3.0 M_⊙ and M_Ac2sin³i = 40.7 ± 2.2 M_⊙.

The integral spectrum of the system was classified O9.7 Ib by Walborn (1972) and Sota et al. (2014). Parkin et al. (2011) studied the Chandra satellite X-ray spectra secured over about 2 years and tried to model them in terms of colliding winds. These authors found some support for the continuing mass exchange between the components Ac1 and Ac2. They also revised the component properties and gave a new spectral classification O9.7 I for component Aa1, and O8 III for component Ac1².

A detailed study of the system was published by Walker et al. (2017). These authors derived new RVs and orbital solutions for both orbits and also new photometry and the light-curve solution. Using these, Walker et al. (2017) estimated improved values of system masses, radii, and effective temperatures. They also concluded that the period of mutual orbit of the two binaries Aa and Ac is close to 49 years.

Mayer et al. (2020) published a conference report on the partial results of the study presented here. They obtained the first RV curve of component Ac2 from the spectra in the optical region and, in combination with photometry, derived the masses of both bodies of the eclipsing subsystem, finding that Ac2 is the more massive of the two. Mayer et al. (2020) also concluded that the system is a semi-detached one. Analysing variations due to light-time effect and secular changes of the systemic velocities of both, Aa and Ac systems, they estimated the orbital period of their mutual revolution to be about 11 700 d (32 years).

Blackford et al. (2020) analysed numerous sets of photometric observations of QZ Car including two series of BRITE satellite photometry secured in 2017 and 2018. They derived epochs of many binary mid-eclipses and constructed a new O − C diagram covering 48 years, using their linear ephemeris

(1)

The authors concluded that the Aa-Ac binary had not yet completed one full orbit since the first set of its photometric observations in 1971 and noted cycle-to-cycle changes in the shape of the orbital light curve with the 5999 period.

An overview of the orbital elements published by various investigators is provided in Appendix A: In Table A.1, solutions for the Aa1–Aa2 207 spectroscopic binary subsystem are listed; those for the eclipsing 60 subsystem Ac1–Ac2 are in Table A.2. Several published estimates of the basic physical properties of all four components are in Table 3.

3. Observational material used in this study

3.1. Spectroscopy

We compiled all published RVs, derived RJDs whenever necessary, and analysed 79 new electronic spectra from three echelle spectrographs. A log of all observations available to us is presented in Table 4. Details of the reduction are provided in Appendix B. There are some problems with two of our principal data sets. Echelle spectra from the CTIO Chiron echelle spectrograph have a very steep blaze function in each echelle order and an accurate spectral normalisation is not easy, especially in cases where strong spectral lines are located near the peak of the echelle order or are close to the end of it. The BESO spectra have low sensitivity, and therefore a rather poor signal-to-noise ratio (S/N), and are useful for our purpose only longward of some 4800 Å. The best set of available spectra is the collection of eight FEROS spectra, all of which have a S/N over 200 in the whole available spectral range. Their distribution over orbital phases of both close orbits is rather good. These are the only spectra that can be investigated in the spectral region from about 4000 to 4400 Å, where one finds useful spectral lines sensitive to the radiative parameters of the system components. To have a good starting point for the comparison with published RVs, we measured RVs in all our new spectra using the program SPEFO (Horn et al. 1996; Škoda 1996), namely the latest version 2.63 developed by J. Krpata (Krpata 2008). It is not easy to find suitable spectral lines for the RV measurements. Blends, the presence of telluric lines in the red parts of spectra, and rather strong interstellar lines all complicate the task. After some trials, we gave up on more sophisticated methods of RV determination. Guided by previous experience, we finally selected two stronger singlet He I lines at 4922 and 5016 Å, and He II 5414 Å. In Fig. 2 we show a selection of the He I 4922 Å and He I 5016 Å line profiles for some spectra, where a weak line of component Ac2 is visible. Whenever possible, we tried to also measure its RV.

Fig. 2. Three examples of the spectra, in which the lines of three components of the system, Aa1, Ac1, and Ac2 can be measured. Top: He I 4922 Å line. Bottom: He I 5016 Å line. RJDs of mid-exposures are shown.

3.2. Photometry

We compiled all available accurate photometric observations and tried to reduce them onto a comparable system of magnitudes that is close to the Johnson photometric system. Details of the reduction are provided in Appendix C and a journal of the data is shown in Table 5. All homogenised ground-based photometric observations are presented in Table 6, which is available in electronic form only.

4. Subsequent data analysis

In spite of an extensive set of observational data, the analysis of a multiple star system as complicated as QZ Car is difficult and one has to proceed in consecutive steps. We first used all homogenised photometric data in the V band and all available BRITE and TESS space photometries to derive the most accurate mean orbital period of system Ac and the epoch of the primary minimum. A formal solution of all these observations in the program PHOEBE1 (Prša & Zwitter 2005, 2006) led to the following photometric mean ephemerides:

(2)

(3)

The phase plot of all V data, TESS, and BRITE photometries for this ephemeris is shown in Fig. 3. The V band data spanning the longest time interval clearly show that there are minor changes in the times of minima, which are related to the light travel time effect in the wide orbit, as we show below. Additionally, there are obvious cycle-to-cycle changes in the light-curve shape, which cannot (only) be due to possible microvariability of the comparison stars used, because they are clearly seen in all space photometries, as shown in Fig. 3. These phenomena inevitably affect the accuracy with which the times of minima and the elements of the light-curve solution are derived.

Fig. 3. Light curve of the eclipsing system Ac based on all homogenised V band observations at our disposal, and on all BRITE R band and TESS broad band space photometries plotted for the ephemeris (2). We deliberately show the BRITE data for their original ‘setups’, which have different zero levels (see Appendix C for details). We note that the bottom two BRITE light curves have lower amplitudes than the rest. This is caused by the fact that part of the stellar light was captured by pixels in the defective part of the detector. These light curves are not usable for the final light-curve solution but can still be used to make local determinations of the times of minima discussed here.

We then proceeded to new orbital solutions based on all available RVs in order to obtain improved ephemerides for both orbits. For all new orbital solutions presented here, we used the program FOTEL (Hadrava 1990, 2004).

4.1. Orbit of the spectroscopic component Aa1

To derive an updated and improved solution for component Aa1 we used all compiled and new SPEFO RVs. Inspection of SPEFO RVs versus phase of the 2073 period led to a surprising finding. While the RVs for the He I 4922 and 5016 Å both define the orbital changes quite well, the RVs from the He II 5411 Å line showed no dependence on the phase of the 2073 period. The mean RV of this line is about +20 km s⁻¹, significantly more positive than the systemic velocity of the whole system. This is clearly seen in the trailed spectra created with a Python program of J. Jonák and shown in the left panels of Fig. 4. We therefore adopted the mean RVs of the two He I lines only to study the 2073 orbital changes.

Fig. 4. Trailed spectra of all FEROS, BESO, and CHIRON line profiles of the three He spectral lines versus orbital phase, for the 2073419 period (left), with phase zero at RV maximum from solution C of Table 7, and for the 5998682 period (right), with phase zero from ephemeris (2).

In order to also monitor the slow change in the systemic velocity as the subsystem Aa moves in a wide orbit with the pair Ac, we sorted the data in time and assigned individual γ velocities to data subsets for time intervals not longer than some 1100 d and much shorter in most cases. Our solution is shown in Table 7 and the corresponding RV curve is in Fig. 5. Individual local γ velocities, together with their mean RJDs, standard errors, number of individual RVs and time interval covered, are provided in Table D.1. For photographic spectra from spectrographs 3 and 4 (see Table 4), we used the published measurements for different groups of lines and treated them as different spectrographs. This confirmed the existence of some significant differences in the RV zero point between different groups of lines, as already noted by previous investigators.

Fig. 5. Radial velocity curve of component Aa1 plotted for the period of 2073419. Top: older data prior to RJD 45000. Bottom: more recent data for RJDs after 48000.

Splitting the data into two subsets and calculating trial orbital solutions for them, we verified that there is currently no compelling evidence of either a secular period change or apsidal motion. However, it turned out that –probably due to a generally lower resolution of older, photographic spectra, and a relatively small RV amplitude of component Aa1– there is a significant difference in the semi-amplitude between the older and more recent electronic spectra. Therefore, we derived two other orbital solutions, for the older and more recent data, keeping the orbital period fixed at the value from the solution for all data. These are also listed in Table 7 and the final local systemic velocities in Table D.1 are based on these two separate solutions. The corresponding RV curves are shown in Fig. 5.

To investigate the character of relatively strong and variable Hα emission, we also produced trailed spectra for all FEROS and CTIO Hα line profiles versus the phases of both the 207 and 5999 periods. These are shown in Fig. 6 where we can see that the bulk of the emission moves in orbit with component Aa1, although it is possible that a part of it originates in other places within the system, probably being associated with the circumstellar matter related to mass transfer from component Ac1 to Ac2. No trace of the secondary Aa2 is seen in the upper panel of the plot.

Fig. 6. Trailed spectra of the Hα line profiles for all FEROS and CTIO spectra plotted versus phase of the 207 orbit with phase zero at RV maximum of the component Aa1 from solution C of Table 7 (top panel), and versus the phase of the 5999 period for ephemeris (2) (bottom).

4.2. Orbit of the eclipsing components Ac1 and Ac2

4.2.1. Orbital solution, mass ratio, and radiative properties

We first used all compiled and new SPEFO RVs of component Ac1 and it was possible to use the mean RV of all three blue lines, He I 4922 and 5016 Å, and He II 5411 Å. We also verified that the RV curves derived from older and more recent data have the same (large) semi-amplitude; therefore a single orbital solution for the determination of local systemic velocities was used. This solution is provided in Table 8 the corresponding RV curve is in Fig. 7, and the local systemic velocities with their rms errors and the range of the spectra are summarised in Table D.2. It is encouraging that this solution led to a mean orbital period, which is consistent within the derived errors with that from ephemeris (2), obtained from fully independent photometric data. Values of the epoch of the primary minimum from both solutions differ slightly more but this is not surprising given the presence of the light-time effect and the different time distribution of photometric and RV data.

Fig. 7. Radial-velocity curve of component Ac1 plotted with the ephemeris of Table 8. We note that component Ac1 is eclipsed in the secondary minimum (phase 0.5).

The less numerous (and much less accurate) RVs of component Ac2 cannot be directly used for local γ determination. However, we note that there are two useful series of spectra, secured within relatively short time intervals, for which we derived solutions with the orbital period fixed and allowing for the mass ratio determination. The resulting mass ratios M_Ac2/M_Ac1 were

FEROS+BESO over RJD 53739–54957: 1.29 ± 0.18

and

CTIO over RJD 58868–58908: 1.29 ± 0.16.

This provides the best estimate of the mass ratio of the eclipsing pair so far and confirms that the brighter, mass-losing component Ac1 is the less massive of the two, contrary to what Stickland & Lloyd (2001) found from the IUE spectra.

We carried out yet another test using eight superb FEROS spectra in the blue spectral region from 3990 to 5440 Å (see one such spectrum in Fig. 8). We used the program PYTERPOL (see Nemravová et al. 2016 for details) to estimate the radiative properties, relative luminosities, and also RVs of components Ac1, Ac2, and Aa1. We note that Sanchez-Bermudez et al. (2017), using different modelling approaches, obtained relatively precise estimates of the flux brightness ratio of the two close binaries F_Ac/F_Aa (see their Table 3). These imply the following estimates of the relative luminosity L₃ of component Aa: 0.585 ± 0.007 for the Markov chain Monte Carlo (MCMC) fit, 0.60 ± 0.04 for the LitPro fit, 0.546 ± 0.001 for the CANDID fit (Gallenne et al. 2015), and 0.57 ± 0.04 for the average fit. Sanchez-Bermudez et al. (2017) adopted the average of the best fits by individual methods. We are not sure that this is the best approach, considering the different levels of sophistication of these models. In our opinion, the CANDID model should be adopted, which also gives the lowest normalised χ². During the calculations, we omitted short segments around Hγ and Hβ lines, and the interval between 5130 and 5350 Å, affected by small emission features. Via artificial rectification, we also removed stronger diffuse interstellar bands near 4428, 4726, 4762, and 4780 Å. Uncertainties of the fit were estimated with the MCMC modelling. The results are in Table 9. The radial-velocity curves of components Ac1 and Ac2 derived by PYTERPOL are shown in Fig. 9. The corresponding orbital solution gives a mass ratio M_Ac2/M_Ac1 = 1.25 ± 0.14 and K_Ac1 = 255 ± 15 km s⁻¹, in good agreement with the results based on all spectra.

Fig. 8. Blue part of the FEROS spectrum taken on RJD 54603.5911, with lines of components Ac1 and Ac2 near one elongation.

Fig. 9. Radial-velocity curve of components Ac1 (blue) and Ac2 (red) based on the PYTERPOL fit of the blue FEROS spectra. Residuals are shown by small dots.

4.2.2. Combined light-curve and RV-curve solution for system Ac

We derived a final PHOEBE solution using all calibrated U​BV data sets, one TESS and one well-covered BRITE light curve, and adopting the effective temperatures of 31 980 K for component Ac1, and 34 300 K for component Ac2 from the PYTERPOL solution (see Table 9). During each iteration, we modified the synchronicity parameter of component Ac2 (the ratio between orbital period and the period of rotation estimated from the instantaneous values of the radius, inclination and v sin i derived with PYTERPOL), the final value being about 3.9. The corresponding solution is shown in Table 10. The system is semi-detached, and the component Ac1 is filling the Roche lobe and is obviously sending material towards Ac2.

The relative contribution of the third light (i.e. flux from the subsystem Aa, which is predominantly from Aa1) is decreasing towards shorter wavelengths, confirming that the effective temperature of component Aa1 is lower than those of components Ac1 and Ac2. There is relatively satisfactory agreement between the results from PYTERPOL and our PHOEBE solution. The fit of the light curves used is shown in Fig. 10.

Fig. 10. Fit of U​BV (top) and BRITE and TESS light curves (bottom) with PHOEBE1. Phases are calculated with ephemeris (2).

The combined solution of the calibrated light curves and a RV curve of component Ac1 for the fixed mass ratio 1.29 presented in Table 10 provides the following standard U​BV values for individual components (if one uses the U​BV values of the whole system in elongations defined by the model fit):

Ac1: V = 7325, B − V = 0029, U − B = − 0969,

Ac2: V = 8615, B − V = 0008, U − B = −0983,

Aa1: V = 6940, B − V = 0251, U − B = −0676.

A standard dereddening for the components of the eclipsing subsystem Ac then gives E(B − V) of 0357 and 0337 for components Ac1 and Ac2, respectively. Using the mean of these two estimates, E(B − V) = 0347, one obtains

Ac1: V₀ = 617, (B − V)₀ = −032, (U − B)₀ = −122,

Ac2: V₀ = 753, (B − V)₀ = −034, (U − B)₀ = −124.

Together with the adopted effective temperatures and bolometric magnitudes and using Flower (1996) bolometric corrections, one obtains M_V = −611 and −485 for Ac1 and Ac2, respectively, which implies a distance towards the system of about 2860–3000 pc, in broad agreement with the distance of the Collinder 228 cluster.

If component Aa1 is a genuine supergiant, then for its estimated T_eff = 28 900 K, it should have a deredenned (B − V)₀ = −025 according to calibration by Flower (1996). For the above-quoted deduced B − V = 0251, this would imply a significantly higher E(B − V) of about 050.

We note that the parallax of QZ Car as given in Gaia DR2 and EDR3 (without account of the multiplicity of the object) seems to be overestimated. The surrounding close-by stars from the cluster have parallaxes that are about half those of QZ Car itself. The distance to the cluster was estimated from 29 O stars by Shull & Danforth (2019) to be of 2.87 ± 0.73 kpc from Gaia DR2, and 2.60 ± 0.28 kpc from photometry, which are similar to what we estimate from the combined solution above. However, it must be kept in mind that nebular emission is known for the Carina nebula and that the observed reddening could be partly affected by that or even by the circumstellar matter within QZ Car itself. Therefore, it may not provide a good estimate of the distance to the system.

The radius of component Ac1 and the orbital inclination from the combined solution presented in Table 10 would imply v sin i = 161 km s⁻¹ for this component, while PYTERPOL gives 126 ± 14 km s⁻¹. The combined solution also predicts log g = 3.27 [cgs] for component Ac1, while the PYTERPOL value is 3.51 ± 0.06. The agreement in this case is not ideal.

4.3. Wide orbit

The current estimates of the mutual orbital period of the two pairs, Aa and Ac, are still rather uncertain, ranging from about 32 to 50 years (11 700 − 18 300 d). There were several attempts to estimate the orbital elements of the long orbit, mainly based on the eclipse timing variation (ETV) analysis (Mayer et al. 2001; Walker et al. 2017; Blackford et al. 2020). However, as only one minimum and one maximum of the O − C diagram have so far been covered (see the top panels of Fig. 11), there is a rather large range of possible values of the period of the mutual orbit of pairs Aa1–Aa2 and Ac1–Ac2.

Fig. 11. A combined fit of data for the long period of 14 403 d. Top panel: The O − C deviations from the times of primary minima. Middle panel: systemic velocities of components Aa1 (black) and Ac1 (red). Bottom panel: corresponding spatial long orbit, which fits the four precise interferometric observations, with the pairs Aa and Ac resolved.

Having the richest and the most homogeneous data sets, we attempted to obtain a more accurate estimate of the elements of the wide orbit. To demonstrate the inherent uncertainties in the determination of the times of binary eclipses, we derived them in two different ways:

(1) We derived local estimates of the times of the primary eclipse of the Ac pair using a formal light-curve solution with PHOEBE1 over limited time intervals and keeping the orbital period fixed at the value from the complete solution; see ephemeris (2). This may provide a more robust estimate of the minima than the estimates for the individual, often not ideally covered minima. The results are provided in Table D.3. We then derived the eclipse timing residuals (ETR) with respect to ephemeris (2).

(2) Alternatively, the individual light curves were also analysed in the following, perhaps more customary way. Individual photometric datasets were used separately as input for the semi-automatic routine automatic fitting procedure (AFP; see Zasche et al. 2014). This routine returns individual times of the primary and secondary eclipses using a light-curve template provided by PHOEBE. It therefore provides up to twice as many eclipse times. However, their accuracy depends on how well the individual minima are covered by observational data and it can often be slightly poorer than the accuracy of the minima defined by the whole local light curves. For that reason, only reasonably well covered minima were used. This is why the corresponding ETRs show a lower scatter than those derived using the first method. The instants of the primary minima derived in this way are shown in Table D.4, and those for the secondary minima are shown in Table D.5.

We attempted to estimate the period and some other elements of the wide orbit for the ETRs derived using both of the above-mentioned ways. We modelled individual quantities formally as a motion in an eccentric orbit using FOTEL. These are listed as solutions 1 and 2 in Table 11. As both these solutions led to similar values of the long orbital period and eccentricity, we combined both sets of minima and derived a joint final solution, which is presented in Table 11 as solution 3.

We then derived orbital solutions for the wide orbit using the local systemic (γ) velocities of both components, Aa1 and Ac1. These are also listed in Table 11 as solutions 4 and 5. A free convergence of all elements for the local γ velocities of component Aa1 led to an unconstrained value of orbital eccentricity of 0.62 ± 7.8, and so we derived solution 3 for Aa1 with eccentricity fixed at the value obtained from both ETRs fits.

Our results illustrate the remaining, non-negligible uncertainties of the current knowledge of the elements of the long orbit but show that the orbital period is probably close to about 14 500 d. Very regrettably, no RVs of components Aa1 and Ac1 are available for the time interval from RJD 51000 to RJD 53000, where there should be the largest difference between γ velocities of systems Aa and Ac. The value of the orbital amplitude of component Aa in the wide orbit thus strongly depends on the early photographic spectra. From the results in Table 11, one can only conclude that the masses of systems Aa and Ac are of the same order of magnitude, system Aa being probably more massive than Ac.

There is also another, and potentially very important source of information, namely the new interferometric observations obtained recently by Sana et al. (2014) and Sanchez-Bermudez et al. (2017), who directly resolved the two pairs. Considering the large uncertainty of the observation obtained with the Sparse Aperture Masking (SAM) instrument by Sana et al. (2014), we decided to skip this one and we used only the two other data points. Moreover, we used two more recent interferometric observations (from March 14, 2017, and April 27, 2017) secured with the GRAVITY instrument attached to Very Large Telescope Interferometer (VLTI), which we found in the ESO archive³. These four positional measurements and our complete sets of the ETRs discussed above were used to get a combined solution for the long orbit. As only four positional observations secured over a limited time-interval are insufficient to define the wide orbit reliably, even more so given that all of them exhibit a non-negligible correlation between the separation and the position-angle error bars, we derived a constrained solution adopting the period and eccentricity from the light-time effect fit, that is solution 3 of Table 11. Corresponding ETRs are compared to the light-time effects and to systemic velocities in the top and middle panels of Fig. 11. The light-time effect due to the mutual orbit of the two pairs is clearly visible, as is the antiphase variation of the γ velocities of components Aa1 and Ac1. The projection of the orbit is shown in the bottom panel of Fig. 11, where the error ellipses are also plotted (we note that for the second one, the observation from GRAVITY (Sanchez-Bermudez et al. 2017), the ellipse is so small that it cannot even be seen around the black dot).

The solution gives an inclination of the wide orbit i₃ = 96° 8 ± 8° 1, Ω₃ = 147° 6 ± 6° 5 and requires a parallax of 000045, implying a distance to the system of over 2200 pc. The mass function of 14.17 then implies the mass of the whole Aa system of 57.1 . Using the observed mass function for the 207 spectroscopic Aa1-Aa2 orbit from solution C of Table 7 and assuming approximate coplanarity of all three orbits, one obtains M_Aa1 = 47.6 and M_Aa2 = 9.5 . Our nominal model thus leads to the total masses of systems Aa and Ac, which are almost equal.

5. Discussion of the emerging picture of the quadruple system

We made attempts to disentangle the eight good-quality blue FEROS spectra using the program KOREL (Hadrava 1995, 1997). We were unable to find a unique and stable solution, which is unsurprising given that those spectra only cover a limited part of the mutual orbit of the Aa and Ac pairs. However, there are indications that the emission components of Hβ and Hγ lines are associated with the mass-gaining component Ac2. This in turn means that what we actually observe for Ac2 could be a pseudophotosphere, that is, inner optical parts of the accretion disk.

We also tried to analyse the residuals from the light-curve solutions, separately for the V magnitude data, and TESS and BRITE photometries. The only detected frequencies were close to the one-year alias of the orbital period and its harmonics. We conclude that the obvious changes in the shapes of the light curves from one cycle to another are stochastic by nature, probably related to changes in the circumstellar matter. Similar changes are known for β Lyr (Larsson-Leander 1970; Harmanec et al. 1996) and some other interacting binaries.

There is also some indication of a slight secular change of the total brightness of the object in the standardised yellow magnitude but confirmation will require future accurate differential observations relative to a non-variable comparison star. For completeness, we also mention that the residuals from the fit of light-time changes based on local fits of complete light curves (solution 2 of Table 11) can be reconciled with a period of 651 d.

Current knowledge of the physical properties of the Aa system is still rather uncertain because we do not see the spectrum of its secondary Aa2. Consequently, there is no RV curve of the companion and no direct estimate of the binary masses. Also, the orbital inclination is not known. We only have the mass function and can estimate the system properties if we adopt some value of mass for component Aa1.

6. Conclusions

We collected and critically evaluated a large body of observational data and obtained much improved orbital solutions for subsystems Ac and Aa of the nine-component system HD 93206. We estimated the mass ratio of the Ac system to be 1.29, confirming that the more massive object is a mass-gaining and less luminous secondary, which is eclipsed in the primary minimum. This system is obviously in the still rather rapid phase of mass transfer between the components, and is observed relatively shortly after the mass-ratio reversal. It is also conceivable that even the Aa1 component of the 207 spectroscopic binary is actually a shell star, a product of previous mass transfer in that system.

Our analysis, in the context of what is currently known, shows that a better understanding of this unusual system will require not only new series of high-S/N spectra covering also the blue part of the optical spectrum but also more sophisticated models. A first step in this direction is the accompanying Paper II (Brož et al. 2022).

Acknowledgments

We thank J.A. Nemravová, who left the astronomical research in the meantime, for her contributions to this study and for the permission to use her Python programs, PYTERPOL and several auxiliary ones. We also profitted from the use of public versions of the program PHOEBE 1 by A. Prša and T. Zwitter, and the programs FOTEL and KOREL, written by P. Hadrava. We acknowledge the use of photometric observations secured by the late Harry Williams at his private Milton Road Observatory and the FRAM V photometry kindly provided by M. Mašek. The research of P.H., M.B., M.W., and J.J. is supported by the grant GA19-01995S of the Czech Science Foundation. B.B. is supported by the NSF grant AST-1812874. J.L.-B. acklowledges support from FAPESP (grant 2017/23731-1). H.B. acknowledges financial support from the Croatian Science Foundation under the project 6212 “Solar and Stellar Variability”. This work has made use of data from the European Space Agency (ESA) mission Gaia (https://www.cosmos.esa.int/gaia), processed by the Gaia Data Processing and Analysis Consortium (DPAC, https://www.cosmos.esa.int/web/gaia/dpac/consortium). Funding for the DPAC has been provided by national institutions, in particular the institutions participating in the Gaia Multilateral Agreement. The paper includes data collected by the TESS mission, which are publicly available from the Mikulski Archive for Space Telescopes (MAST). Funding for the TESS mission is provided by NASA’s Science Mission directorate. We gratefully acknowledge the use of the electronic databases: SIMBAD at CDS, Strasbourg and NASA/ADS, USA. We also acknowledge the use of the WDS catalogue of visual binaries and multiple systems and thank to Brian Mason for providing us with individual astrometric observations of the components of HD 93206 that are collected in the catalogue.

References

Appendix A: Published orbital solutions for the two subsystems

The orbital solutions for the spectroscopic 207 Aa1-Aa2 binary are summarised in Table A.1, and those for the eclipsing 60 Ac1-Ac2 binary are shown in Table A.2.

Appendix B: Details of the spectral data sets, reduction and measurements, and comments on archival RVs

Initial reductions of FEROS and BESO spectra were carried out by A.N., and those of CTIO CHIRON spectra were carried out by R.C.-H. and W.F. Normalisation, removal of residual cosmics, and flaws and RV measurements of all spectra were carried out within the program SPEFO (Horn et al. 1996; Škoda 1996), namely the latest version 2.63 developed by J. Krpata (Krpata 2008). SPEFO displays direct and flipped traces of the line profiles superimposed on the computer screen that the user can slide to achieve a precise overlapping of the parts of the profile to be measured. Following our usual approach, we also measured a selection of good telluric lines in the red spectral regions to obtain an additional fine correction of the RV zero-point of each spectrogram. The wavelength zero-point of BESO spectra was also checked via Na I 5889 and 5895 Å interstellar lines.

Here we provide comments on the RV data sets compiled from the astronomical literature, which are listed in Table 4. File A: These RVs are from the Radcliffe 2-prism photographic spectra secured in 1956, calibrated to the IAU RV standards. All RVs of QZ Car are denoted as uncertain by the authors. Two of the four spectra were measured repeatedly by different measurers, revealing RV differences of up to 14 km s⁻¹. It is not a priori clear to which component they refer. They might be reconciled roughly with the 20734 orbit of component Aa1, but with a large scatter and very negative mean RV. Plotting them versus phase of the 59987 period resulted in a scatter diagram. The Ca II interstellar lines were also measured, their mean RV being +0.4 ± 1.6 km s⁻¹ from six measurements.

File B: These three photographic RVs have a small RV range from −11 to 0 km s⁻¹ and do not follow either the 20735 or 59987 orbit, as already demonstrated by Mayer et al. (2001); see their Fig. 5. Their mean RV near RJD 37557 is −2.3 km s⁻¹. As a precaution, we do not use them in our analyses.

File C: These near-UV/blue and yellow/red photographic spectra from Cerro Tololo Observatory show the differences in the semi-amplitudes and systemic velocities of RV curves of different spectral lines. For a better comparison with our new spectra, we used the He I RVs, but even so one has to be aware that the RV zero-point might differ from other sources due to line averaging and blending effects. The mean RV of the interstellar Ca II lines for these spectra is −1.6 ± 1.2 km s⁻¹. Although the tabulated JDs are denoted as heliocentric ones by the authors, they are regrettably only given to two decimal places (we note that the JDs and HJDs start to differ only from the third decimal place onwards).

File D: These photographic spectra were secured at the coudé focus of the ESO 1.5 m telescope in 1977-1978 and have a linear dispersion of 12 Å mm⁻¹. RVs were measured using the ARCTURUS oscilloscopic measuring machine but the tabulated RVs are the mean RVs of components Aa1 and Ac1 for He I 3820, 4026, and 4472 Å, and Si IV 4089 Å lines. Mean RVs of Hγ, H8, H9, H10, and H11 are also tabulated. The rms errors of the mean RVs are rather high, which is probably due to the effects of line blending and differences in the RV amplitudes. The RV zero-point of these mean RVs is therefore rather uncertain. They do not cover the minimum RV of the 207 orbit. Also in this case, Julian dates denoted as heliocentric ones are given to only two decimal places.

File E: These RVs are nine cross-correlation RVs from the SWP high-resolution IUE images over the range 1250 to 1900 Å. Their zero-point of the RV scale was tuned via a mean interstellar velocity from three sources, and therefore they seem usable for the monitoring of the systemic velocity changes.

File F: These RVs come from the first available electronic spectra secured with the spectrograph Echelec attached to the ESO 1.5 m reflector in 1992 and 1993, and the coudé spectrograph CES attached to the coudé auxiliary telescope (CAT) of the 1.4 m ESO reflector in 1994. The Echelec spectra mainly cover the Hβ and He I 4922 Å lines. The CAT spectra contain only the neighbourhood of the He I 4922 Å line and two of them were centred on Hα. The authors warn that the Hβ line is broad and less suitable for RV measurements and that the 1992 spectra are affected by a flaw between 4914 and 4917 Å. Regrettably, we were unable to obtain the original spectra and the RV zero-point is somewhat uncertain.

File H: These RVs come from the high-resolution Cantenbury University HERCULES spectrograph and are based on a selection of symmetric, mainly He I lines but also some high-excitation lines.

Appendix C: Details of the photometric data files and their homogenisation

We did our best to reduce the yellow, blue, and ultraviolet photometric observations from various sources into a system that is similar to the standard Johnson system. This is not an easy task because the observations available to us consist of a mixture of differential observations measured relative to various comparison and check stars, which are listed in Table C.1, and all-sky observations.

As the three of the most frequently used comparisons were observed with the TESS mission (Ricker et al. 2015), it is possible to provide some comments on their photometric stability: HD 93131 = TIC 390670978: The TESS signals are mainly that of stochastic low-frequency variability. The object is a Wolf Rayet (WR) star (WN6ha-w; Hamann et al. 2006). WR stars are expected to exhibit stochastic low-frequency signals (likely due to wind clumps), a good relatively large-amplitude example being WR40 studied by Ramiaramanantsoa et al. (2019). The strongest signal in the frequency spectrum of HR 93131 (which does not stand out above the red noise ‘forest’ of the low-frequency stochastic variability) is at f = 0.365697 c d⁻¹, with an amplitude of 3.25 ppt. However, the maximum–minimum (max–min) amplitude is 3.5% (ten times higher than the strongest peak in the frequency spectrum). As this variability is stochastic, it is not expected to introduce any artificial signals in the target star when differential photometry is performed. The main effect of this genuine variability seen in TESS is an effective increase in the noise of the ground-based measurements made relative to HD 93131. The standard deviation of the relative flux measurements of TESS for this target is 0.0067.

HD 93695 = TIC 459164314: There are multiple low-amplitude periodic signals in the TESS data for this system, plus some low-amplitude low-frequency stochastic excess. The strongest five signals are f₁ = 0.123124 c d⁻¹ (amplitude = 0.83 ppt), f₂ = 0.247167 c d⁻¹ (clearly 2f₁, amplitude = 0.77 ppt), f₃ = 0.56325 c d⁻¹ (amplitude = 0.65 ppt), f₄ = 1.291885 c d⁻¹ (amplitude = 0.2 ppt), and f₅ = 4.829406 c d⁻¹ (amplitude = 0.16 ppt). The max–min amplitude of the TESS relative flux is 0.57%, and the standard deviation is 0.0011. Published spectral classifications of HD 93695 range from B3 (Houk & Cowley 1975) over B5 (Feast et al. 1955) and B8/9 (Loden et al. 1976). Feast et al. (1955) also reported rather large RV variations for this star. Figure C.1 shows a phase plot of the RVs of these latter authors versus the phase of the 812189 period, which corresponds to TESS frequency f₁. It is obvious that HD 93695 is a new spectroscopic binary, and probably has a slightly eccentric orbit. As a free solution diverges to an unrealistically high eccentricity because of an incomplete phase coverage, we derived only a trial circular orbit solution, which gives K = 65.6 ± 8.0 km s⁻¹, an epoch of maximum RV of RJD 34095.37 ± 0.15, and a systemic velocity of +11.6 ± 6.3 km s⁻¹. The first CHIRON CTIO echelle spectrum of the star was obtained on November 26, 2021, and it clearly shows the presence of a sharp-lined B-type secondary, which is seen in all He I lines. A segment of this spectrum in the neighbourhood of the He I 6678 Å line is shown in the bottom panel of Fig. C.1. A detailed study of this binary will be the subject of a separate paper. Obviously, the differential observations reduced relative to this comparison should be treated with caution. They fortunately represent only a small fraction of the data we are using here.

Fig. C.1. Top: Phase plot of Feast et al. (1955) RVs of HD 93695 for the period 8121289, which corresponds to the strongest f1 frequency we find in the TESS photometry. The epoch of phase zero is RJD 34095.37. Bottom: Segment of one CHIRON spectrum of HD 93695 taken on November 26, 2021, on which all He I lines are split into one broad and one sharp component. This proves that the object is SB2 composed of two B stars.

HD 93222 = TIC 391050404: As for HD 93131, the TESS variability appears to be purely that of low-frequency stochastic excess. The strongest signal in the frequency spectrum is at f = 0.1672283 c d⁻¹ (amplitude = 0.8 ppt). The max–min amplitude of the normalised TESS flux is 0.87%, and the standard deviation is 0.0014. Such low-amplitude stochastic variability is typical for O stars viewed with TESS (Bowman et al. 2020). The star is listed as having a spectral type of O7V((f))z by Sota et al. (2014).

Here, we provide details on the individual data sets identified by their station numbers from Table 5. Station 11: These U​BV observations were obtained in 2004 by Marek Wolf and in 2011 by Petr Zasche. They were transformed to the standard Johnson system via non-linear transformations with the latest release of the program HEC22 that is able to model time variability of the extinction during the nights (Harmanec et al. 1994). Station 12: These U​BV observations were obtained by Pavel Mayer, Horst Drechsel, and Reinald Lorentz and were later also re-reduced with the program HEC22 to the standard system (Mayer et al. 1992, 2010). Station 61: These Hipparcos H_p observations were transformed to Johnson V magnitude after Harmanec (1998) assuming B − V = −0140, and U − B = −0832 for QZ Car. This seems justified because the B − V and U − B colours of QZ Car show almost no variability with either the orbital phase of the Ac system or with time. Only data with error flags 0 and 1 were used. Station 93: These all-sky ASAS3 Johnson V observations were processed in several steps. First all observations of grade D were omitted. As the data are published for five different diaphragms, we selected the one with the lowest average rms error, diaphragh 4 in our case. We then inspected residuals with respect to the light curve and removed all deviating data points. Station 100: These observations were obtained by Mark Blackford at his backyard observatory in Chester Hill, Sydney, between 2015 and 2018. A 0.080 m f/6 refractor and Canon 600D DSLR camera were used. Typically, several hundred images were obtained each night and 5 images were always averaged. The exposure times were 20 s in 2015, 30 s in 2016 and 2017, and 50 s in 2018. These data were re-reduced by M.B. in December 2020 relative to HD 93131, with HD 93695 serving as the check star. Station 105: These observations were secured by Anthony Moffat (Moffat 1977) through a λ5170 Å filter, 190 Å wide. We added the V = 8101 magnitude of his comparison HD 93222 derived from all-sky standard observations at Sutherland by Petr Zasche to obtain approximate V magnitudes of QZ Car. We assigned all these observations with the rms error of 0007 estimated by the author. Moreover, we calculated the heliocentric Julian dates from the published JDs. Station 106: These discovery U​BV observations were published by Walker & Marino (1972). Both the JDs and magnitudes were given to only two decimal places there. A more extended data set was now re-reduced by W.S.G. Walker, corrected for extinction and transformed to the Johnson system. Data between RJD 41290 and RJD 43270 come from observations originally devoted to η Car monitoring. The data were observed in a sequence of three consecutive 10 s integrations, subsequently in V, B, and U, with a 2-5 s delay between integrations in each filter. Julian dates were recorded for the beginning of observations; therefore, to obtain mid-exposure times of these observations, we increased Julian dates by 30+3+15 seconds, i.e. 0.00055 d. For observations with V and B filters only, the respective correction was 0.00037 d. We then derived heliocentric Julian dates for them. We also subtracted the U​BV values of the comparison star HD 93131 derived by Walker in 2019, V = 6496, B − V = −0042, and U − B = −0882 and added the all-sky values of HD 93131 derived at La Silla instead (see Table C.1). Additional U​BV observations with the same instrument were obtained by Grant Christie. A part of these observations was used by Mayer et al. (1998), who incorrectly attributed them to the Mt John Observatory. They were transformed to the standard Johnson system and reduced relative to HD 93695, for which V = 649, B − V = −011, and U − B = −062 was originally used. We re-reduced these observations relative to the standard values of HD 93695 derived by Mayer et al. (2010) at LaSilla (cf. again Table C.1). Station 107: The Bright Star Monitor Australia Station owned and operated by Peter Nelson. This BV​RI photometry was reduced by AH relative to a set of comparison stars in the field of observations. Station 108: Martin Bruce Berry Bright Star Monitor BV​RI photometry, also reduced by AH relative to a set of comparison stars in the field of observations. Station 109: These BV observations were secured by an AAVSO observer Terrence Bohlsen (AAVSO observer code BHQ) in Mirranook Observatory and transformed to the standard Johnson system by him. They consist of three subsets. The observations in 2009 (RJDs 54945.05 - 55043.93) were secured with an SBIG ST9 camera. Later observations were secured with an SBIG ST10XME camera. The observations from 2010 and two from January 2013 (RJD 55218.99 - 56324.01) were reduced differentially also relative to HD 93131 while the 2013 to 2017 observations (RJD 56324.01 - 57837.94) were reduced relative to an ensemble of three comparison stars, HD 93131, HD 93222, and HD 305536 and were reported to AAVSO. Station 110: These observations were secured by the TESS space photometer (Ricker et al. 2015) in sectors 10 and 11 of the mission (between March 26, 2019, and May 21, 2019), corresponding to 54.4 d of near-continuous observation. Light curves were extracted from the Full Frame Images (FFIs) at a 30-minute cadence for use in this work, where the sky-background level was computed and subtracted from each frame by JLB. This detrending method (sky-background subtraction) was preferred over other detrending algorithms which can distort relatively high amplitude and long-timescale variability (e.g. the 6 d binary orbital signal associated with component Ac). Light curves for these sectors are also available with two-minute cadence, which were used to search for periodic signals out to the Nyquist limit of 360 c d⁻¹ (after pre-whitening against the 6 d orbital period and its many harmonics). Besides the 6 d orbital period and its harmonics, no further periodic signals were found (although low-frequency stochastic variability is present). A neighbouring source, HD 305522, is 66″ away from QZ Car and is approximately two magnitudes fainter than it is in the TESS band. We estimate that approximately 4% of the flux in the adopted aperture comes from this neighbour. However, a pixel-by-pixel blending analysis demonstrates that the signals in the extracted light curve can be confidently attributed to QZ Car and that the net effect of blending from HD 305522 (and perhaps other much fainter neighbours) is a very mild suppression of the variability amplitudes. Station 111: These observations were originally obtained as U​BV observations by the late Harry Williams in his private Milton Road Observatory and reduced by W.S.G. Walker. However, for technical reasons (e.g. patina, which affected one of the filter glasses) we could only use the V band observations. Stations 122 and 123: These R band observations were secured with two BRITE nanosatellites (see Weiss et al. 2014), BRITE-Heweliusz (BHr, station 122) in 2017, 2018, and 2020, and BRITE-Toronto (BTr, station 123) in 2018 and 2019. The 2017 and 2018 observations have already been used by Blackford et al. (2020). All data have now been reduced with the latest version of the decorrelation software by A.P. Observations in individual seasons are split into ‘setups’ (see Popowicz et al. 2017 for explanation of what setup refers to here). Observations in a given setup are reduced independently, which means that each setup has a different zero-point in photometry, and magnitude shifts between setups are allowed. However, there are exceptions. If the position of a subraster was the same for different setups of the same satellite, two or more setups were combined, applying magnitude offsets to minimise differences between them. This is the case for station 122 data secured between RJD 57784 and 57806, and station 123 data secured between RJD 58165 and 58187. To reduce the scatter and remove possible systematic errors, we created normal points averaged over the satellites’ orbital periods of 0068194 and 0067431 for stations 122 and 123, respectively. The number of these normal points is quoted in Table 5. Station 124: These CCD V band observations were obtained by Dave Blane with a 150mm f/5 refractor and a Canon 1300D DSLR camera mounted on a GEM goto mounting. Ensemble aperture photometry was facilitated using the IRIS software package. The camera g magnitudes were linearly transformed to the Johnson V band using an MS-Excel spreadsheet. Differential extinction was not applied. Each adopted measurement was the average obtained from ten separate field images. The systematic difference in the mean magnitudes for the two seasons can probably be explained by the fact that different comparison star sets were used. For the first season, the AAVSO values for six comparison stars were used while for the second season, magnitude values for four comparison stars recommended by Mark Blackford were used. In the final step, the V magnitude differences QZ Car − HD 93695 were added to our standard V value for HD 93695 from Table C.1. Station 125: These observations were obtained by Mark Blackford at Congarinni Observatory between 2019 and 2020. Slightly different instrumentation was used in different seasons. In 2020, an 80mm f/6 refractor and Atik One 6.0 CCD camera with Astrodon Johnson V filter were used to record images of QZ Car. Typically, about 40 images with two-second integration times were secured each night (although the number and/or exposure time was different on a few nights). Instrumental V magnitude differences QZ Car minus HD 93131 were derived, HD 93695 serving as a check star. In 2019, the same instrumentation as in 2020 was used but with a 2X tele extender. This meant that HD 93695 was outside the field of view, and so HD 305523 served as a check star instead. Several hundred images were usually obtained each observing night, with 10 second exposure times. Station 126: These instrumental V magnitude observations were obtained by Martin Mašek relative to the comparison star HD 92741 with the FRAM 0.30 m automatic monitor at Pierre Auger Observatory in Argentina (The Pierre Auger Collaboration et al. 2021) and kindly put at our disposal by their author. We added V = 723 to the magnitude differences var-comp.

Appendix D: Supplementary tables

All Tables

All Figures

	Fig. 1. Sketch of the multiple system HD 93206 drawn in arcsec scale on both axes. We note the expanded scale on the ordinate. The red circle is centred on the position of the unresolved Aa1-Aa2 spectroscopic binary. North is up, east is to the left.
In the text

	Fig. 2. Three examples of the spectra, in which the lines of three components of the system, Aa1, Ac1, and Ac2 can be measured. Top: He I 4922 Å line. Bottom: He I 5016 Å line. RJDs of mid-exposures are shown.
In the text

Fig. 3. Light curve of the eclipsing system Ac based on all homogenised V band observations at our disposal, and on all BRITE R band and TESS broad band space photometries plotted for the ephemeris (2). We deliberately show the BRITE data for their original ‘setups’, which have different zero levels (see Appendix C for details). We note that the bottom two BRITE light curves have lower amplitudes than the rest. This is caused by the fact that part of the stellar light was captured by pixels in the defective part of the detector. These light curves are not usable for the final light-curve solution but can still be used to make local determinations of the times of minima discussed here.

In the text

	Fig. 4. Trailed spectra of all FEROS, BESO, and CHIRON line profiles of the three He spectral lines versus orbital phase, for the 2073419 period (left), with phase zero at RV maximum from solution C of Table 7, and for the 5998682 period (right), with phase zero from ephemeris (2).
In the text

	Fig. 5. Radial velocity curve of component Aa1 plotted for the period of 2073419. Top: older data prior to RJD 45000. Bottom: more recent data for RJDs after 48000.
In the text

	Fig. 6. Trailed spectra of the Hα line profiles for all FEROS and CTIO spectra plotted versus phase of the 207 orbit with phase zero at RV maximum of the component Aa1 from solution C of Table 7 (top panel), and versus the phase of the 5999 period for ephemeris (2) (bottom).
In the text

	Fig. 7. Radial-velocity curve of component Ac1 plotted with the ephemeris of Table 8. We note that component Ac1 is eclipsed in the secondary minimum (phase 0.5).
In the text

	Fig. 8. Blue part of the FEROS spectrum taken on RJD 54603.5911, with lines of components Ac1 and Ac2 near one elongation.
In the text

	Fig. 9. Radial-velocity curve of components Ac1 (blue) and Ac2 (red) based on the PYTERPOL fit of the blue FEROS spectra. Residuals are shown by small dots.
In the text

	Fig. 10. Fit of U​BV (top) and BRITE and TESS light curves (bottom) with PHOEBE1. Phases are calculated with ephemeris (2).
In the text

	Fig. 11. A combined fit of data for the long period of 14 403 d. Top panel: The O − C deviations from the times of primary minima. Middle panel: systemic velocities of components Aa1 (black) and Ac1 (red). Bottom panel: corresponding spatial long orbit, which fits the four precise interferometric observations, with the pairs Aa and Ac resolved.
In the text

Fig. C.1. Top: Phase plot of Feast et al. (1955) RVs of HD 93695 for the period 8121289, which corresponds to the strongest f1 frequency we find in the TESS photometry. The epoch of phase zero is RJD 34095.37. Bottom: Segment of one CHIRON spectrum of HD 93695 taken on November 26, 2021, on which all He I lines are split into one broad and one sharp component. This proves that the object is SB2 composed of two B stars.

In the text