Comprehensive Photometric Investigation of an Active Early K-type Contact System -- IL Cnc

Comprehensive photometric investigation was carried out to the early K-type contact binary - IL Cnc. A few light curves from both ground-based telescopes and the Kepler space telescope were obtained (or downloaded) and then analyzed in detail. They are mostly found to be asymmetric and there are even continuously changing O'Connell effect in the light curves from Kepler K2 data, suggesting the system to be highly active. Using the Wilson-Devinney code (version 2013), photometric solutions were derived and then compared. It is found that the calculation of the mass ratio is easily affected by the spot settings. Combining the radial velocities determined from LAMOST median resolution spectral data, the mass ratio of the binary components is found to be $M_2/M_1 = 1.76\pm 0.05$. The components are in shallow contact ($f\sim 9\%$) and have a temperature difference about $T_2 - T_1 = -280\pm 20$ K. The system is demonstrated to be W-subtype, which may be a common feature of K-type contact binaries. The masses of the binary components were estimated to be $M_1\sim 0.51$ M$_{sun}$ and $M_2\sim 0.90$ M$_{sun}$. The values are in good agreement with that deduced from the parallax data of Gaia. The results suggest that the primary component lacks luminosity compared with the zero main sequence. The H$\alpha$ spectral line of the primary component is found to be peculiar. Combining newly determined minimum light times with those collected from literature, the orbital period of IL Cnc is studied. It is found that the (O$-$C)s of the primary minima show sinusoidal variation while the secondary do not. The oscillation is more likely to be caused by the starspot activities. Yet this assumption needs more data to support.


Introduction
K-type contact binaries are important objects because of their special properties (Bradstreet 1985) such as O'Connell effect (unequal heights of light maximum, see O'Connell 1951aO'Connell , 1951b, period cut-off (Paczyński et al. 2006;Rucinski& Pribulla 2008;Li et al. 2019), W-subtype phenomenon (referring to W-and A-type classifications of contact binaries, see Binnendijk (1970)) and etc.. A recent statistical study shows that the period distribution of EW type eclipsing binaries peaks at 0.29 days (Qian et al. 2017) which is quite close to the period range of early K-type contact systems. Combining their special properties, the early K-type contact systems are key to the study of short period contact binaries. With fast rotation and temperature similar to the sun, the early K-type contact systems are supposed to be active. IL Cnc is such an early K-type system which might have these special properties that come into sight. IL Cnc (cross id: GSC 1400-0455, EPIC 211988016, and etc.. α 2000.0 = 08 h 55 m 51 s .49, δ 2000.0 = +20 • 03 38 .6) was first discovered as an EW-type variable by Rinner et al. (2003) based on unfiltered ccd data (Alton 2018). Photometric data were also collected by surveys such as ROSTE-I survey (NSVS, Woźniak et al. 2004) and ASAS survey (Pojmański et al. 2005). With a short period about only 0.267 days and a large amplitude of 0.6 mag, this object was monitored a number of times for minimum light times. Just recently, this target was studied by Alton (2018) through analyzing light curves from the 0.28-m telescope in UnderOak Observatory and the 0.4-m telescope in Desert Bloom Observatory. Using PHOEBE 0.31a, the author found a preliminary result of roughly determined mass ratio about 1.5 to 2.0, but not yet conclusive determined photometric elements, nor any "meaningful" changes in its orbital period. The author also found this target an active system with asymmetric maxima in its light curves. More information about this target has been published thanks to large surveys including Kepler mission (Borucki et al. 2010), LAMOST (Large Sky Area Multi-Object Fiber Spectroscopic Telescope, e.g. Cui et al. 2012;Zhao et al. 2012) and Gaia (Gaia Collaboration et al. 2016). Therefore, it is possible to study this target in more detail.

Photometric Observations from ground-based telescopes
Complete CCD photometric observations of IL Cnc were carried out first on April 04 and 05, 2015, using the 2.4-m telescope of TNO (Thai National Observatory), which is located in Chiangmai, Thailand, and operated by NARIT (National Astronomical Research Institute of Thailand). Four filters B,V,R and I in Johnson Cousin system were used together with the ARC 4K CCD which has a field of view 8.8 × 8.8 and the integration time was 10 s for each image.
In order to make a comparative analysis, the target was subsequently observed on November 16 and 17, 2016, using the sino-Thai 70 cm telescope in Lijiang, China. This time only three filters (VRI) were put to use with a larger field of view about 20.5 × 20.5 , since the Andor technology 2K CCD was used. The exposure times were changed to 60 s (V), 40 s (R) and 30 s (I). Both sets of data were reduced in a standard way using ccdred packages from IRAF. A differential photometry was applied in order to get the net variation (light curves) of the target.
The standard deviations of the C-Ch data were calculated to be 0.021 mag (B), 0.016 mag (V), 0.013 mag (R) and 0.012 mag (I) for the first set, and 0.009 mag (V), 0.007 mag (R),and 0.012 mag (I) for the second set. The original data of these two sets of light curves can be retrieved from online material. The two sets of light curves were phased using the ephemerides min.I (HJD) = 2457118.08544 (8)   Aside from those light curves, a few minimum light times were determined by using a least-squares parabolic fitting method. A few of them were determined by using the data from 60 cm telescope of YNOs in Kunming, China. The newly determined minimum light times are listed in Table 2.

Orbital Period Investigation
This target has been monitored by several investigators over a time span about 20 years. A few minimum light times were published. We collected all the available times of minimum light from literature and database. They are listed in Table 3, including four reprocessed data. These four data were recalculated from ROTSE and ASAS data by using an average method (see Liu et al. 2015a) because they were observed with long cadence and thus sparsely sampled. The  Reference: (1) ROTSE (Gettel et al. 2006 vizier database 3 , a spectral type of K0-K2 was estimated (Cox 2000). Thus, the temperature of star 1 (the component eclipsed at min I) was set to be T 1 = 5000 K (T 1 is usually a little bit higher than the average). This temperature indicates convective envelope for the components.
Accordingly, the gravity-darkening coefficients were set g 1 = g 2 = 0.32 (Lucy 1967) and the bolometric albedo A 1 = A 2 = 0.5 (Ruciński 1969). The square-root functions (LD = −3) were chosen for the treatment of limb-darkening. The corresponding coefficients were calculated by the code according to Van Hamme's table (1993). It is to mention that low-resolution spectra of this target (one example is shown in Figure 3) were obtained by the team of LAMOST, and a K3 spectral type was suggested. Gaia catalogue 4 also lists the temperature of IL Cnc to be 4894 K. Thus, the temperature estimation should be plausible.
Mode 3 (contact model) was assumed initially for this object according to its EW-type light curves and thus the dimensionless potential of components, Ω 1 and Ω 2 were bound to be the same. For other configurations (e.g. Mode 2), the potentials should be set according to the Roche potential limits. Some adjustable parameters in the beginning of the calculation were: the orbital inclination i; the mean temperature of star 2, T 2 ; the monochromatic luminosity of star 1, L 1X (X = B,V,R,I) and the dimensionless potential. The q-search method was applied and solutions were derived for a series of mass ratios. The mean sum of weighted square deviations ( , hereafter mean residuals) along with mass ratios are plotted in Figure 4. For the two sets of light curves, the calculations were started from the same initial parameter values and carried out independently in order to make a comparison.
The minimum of mean residuals were achieved at different values of mass ratio for the two sets (q = 1.5 and q = 2.3 respectively) after the q-search. Then, the mass ratios were set free to start calculation for comprehensive solutions. Here comprehensive solutions refer to that starspots, third lights or other assumptions were taken into consideration. These assumptions are widely accepted in the study of late type short period eclipsing binaries. It should be noted that we don't have bias as to whether to add cool or hot spots to the surface of either component in models because we think any situation is possible unless the data demonstrate it.
The process of calculation was complicated and time consuming while we brought in a grid-search method (see Chapter 8 in the book of Bevington & Robinson (2003)) in seeking for the best solutions. The final solutions are listed in Table 4 for each main case. Parameters for spots are θ s (latitude), ψ s (longitude), r s (radius) and T s /T * (temperature factor), as well as which component they locate on. The parameters which we attempted griding values are marked "trial". "LCs 2015" and "LCs 2016" in the table denote the light curve data from Chiangmai 2.4-m telescope and Lijiang 70cm telescope respectively. It should be noted that for LCs 2015, we tried an evolutionary starspot, which is a new function of W-D code 2013 as mentioned above. The reason we attempted the evolutionary starspot is that we found the light curves have slightly changed during the two continuous nights. This kind of solutions is marked "e-spot" in Table 4. We assume that "e-spot" was seen in the system for only one night and then found the best solutions with "e-spot" present in the first night. This assumption may be somewhat arbitrary, so this set of solutions might be a tentative result. The derived light curves of the best solutions with normal spots (column 2 and 5) and "e-spot" are shown in Figures 5 and 6, respectively. The χ 2 v of each set of solutions are also listed. In order to calculate this parameter, the errors of the light curve data are needed. Combining the errors of C−Ch data, those errors were estimated to be 0.012 mag (B), 0.011 mag (V), 0.009 mag (R) and 0.009 mag (I) for LCs 2015, and 0.006 mag (V), 0.005 mag (R),and 0.009 mag (I) for LCs 2016. The χ 2 v and Σ values imply that for both sets of light curves, solutions with spots are significantly better than solutions without them. It should be mentioned that we tried adding a third light in our calculation but no reasonable third lights were achieved, which suggests the third light should be less than 2-3% of the total luminosity, taking account of the errors.

Anlysis of Kepler data
The Kepler spacecraft and telescope were designed on the main purpose of searching for exoplanets, which however has provided a large amount of high precise and long-term monitored photometric data for eclipsing binaries. IL Cnc was one of them that benefited from this mission in its rejuvenated stage -K2 (Howell et al. 2014 Table 4 were extracted then. The light curves (only a part of lc05 as an example) vs BJD (barycentric Julian Date) and the phased light curves (rebinned, see Section 4.2) are shown in Figure 7.
The light curves used for analysis are PDCsap data. We did not perform further detrending procedure because the general trend of them are flat and detrending may cause further problems.

Information from maximum and minimum light
From Section 3.2, it is known that O'Connell effects present in the light curves of this target. To our expectation, this effect was also found prominent in K2 data. It was analyzed by calculating the difference of two maxima of the light curves in each cycle, which is shown in Figure 8. The maxima were determined with a least squares parabolic fitting method to the combined data (explained in the following paragraph). It is shown clearly that the O'Connell effect changes with time. The shape of the light curve might be changing all the time, which demonstrates that the target is highly active, and so that the "e-spot" scenario in Section 3.2 is plausible.
All the minimum light times determined from K2 data are listed in Table A1. Since there are much less from enough data in each cycle of the light curves, we tried to use the following combining method to determine minimum light times. First, combine data in 10 continuous cycles and fold them in period to determine one minimum (maximum) data, and then shift   13 3 cycles to determine another one. When meeting the interrupted cycles (or the starting and ending parts), make sure there are at least 7 cycles for one measurement. This method helps to ensure there are enough data to determine one data point and avoid large spanning of time.
The O−C data are calculated using the same formula as used in Section 3.1. It should be noted that since the timings of K2 data are in BJD, other data are then unified into BJD 5 as well for the following analysis. The O−C diagram for K2 data is shown in Figure 9. The overall O−C diagram was replotted and shown in Figure 10. Seen from Figure 9, both primary and secondary minima change with time, but they vary differently. Comparing Figure 9 with Figure 8, it is found that the variation of secondary minima in Figure 9 and the O'Connell effect in Figure 8 look quite similar. It might be the case that the times of minimum light are shifted by the spot activities as pointed out by many investigators (Kalimeris et al. 2002;Watson & Dhillon 2004). However, because the primary eclipse is deeper than the secondary (see Figures 1 and 7), the primary minima may not be so seriously affected by the spot activities, as shown in Figure 9. Therefore, we tried to analyze the primary minima separately. As shown in Figure 11, a cyclic variation might be presenting.
The resulted curve does fit the data well (expect for a few data with large errors). The reduced chi-squares of the fit was calculated to be χ 2 v = 1.35. The standard deviation of the O-C residuals is only about 0.00020 days, which is significantly smaller than the amplitude (∼ 0.0013 days) of the cyclic component, which may indicate that the cyclic variation is true for the primary minima.

Photometric solution of Kepler K2 data
The continuous light curve is pretty useful for studying the photometric properties of the target.
However, since the exposure time is 30 minutes, which is rather long compared with its period ∼ 0.267 days, it is necessary to merge data in neighboring cycles. The data were phased and superimposed using Formula 1. Then the phased light curve were rebinned and fitted to derive the new appropriate light curves for synthetic analysis. The magnitude of one point is an interpolation of the phased data spanning about 30 minutes (approximately equals to the exposure time) and there are altogether 200 data points in one light curve. It should be noted that because the shape of lc05 changes greatly (the O'Connell effect is almost reverted, see For the long cadence data, there must be problem of average effect for the data points, which will make sharp features disappear and change the depth of the eclipse. This is the so called smear effect. The W-D code uses the parameter "NGA" to account for this issue and it was set to NGA = 2 in our calculation. For other parameters, similar procedures as mentioned in Section 3.2 were conducted. The best photometric solutions were achieved for each light curve. The results are shown in Table 5. Since both hot and cool spot scenarios led to good fit for lc05a, both results are listed. While for lc05b and lc18, only cool spot scenarios could give good results (Σ of hot spot scenarios are much worse) and hence listed. It is inferred from this table that the spot activities do affect the solutions of the light curves and cause large uncertainties of the mass ratio. However, other parameters are yet consistent with each other, especially the inclinations, which implies the results to be reliable on the whole.

Analysis of Lamost msp data
Except for the Kepler data, the target was also found to be in the field of Lamost survey. This project uses a telescope with effective aperture of 4 meters and a total of 4000 fibers (Liu et al. 2015b), which makes it a powerful tool for spectral acquisition. Our target was observed by the telescope with both low and median resolution modes. Here, we use median resolution spectra (Wang et al. 2019) to study for more detailed information, e.g. the radial velocities (RVs). The detailed information about the spectrograph and the survey could be found in the literature (Cui et al. 2012;Zhao et al. 2012).
There are six coadding spectra of median resolution that were obtained by Lamost group.
However, only three of them are good enough to be used (others have much lower SN ratios).
They (hereafter "msp091","msp119" and "msp151") were observed on nights with local MJD 58091, 58119, 58151, with 10 minutes for each single exposure. Each spectrum is a coadding of three spectra observed at almost the same time. They cover two wavelength ranges: 495-535 nm and 630-680 nm which are the so-called blue band and red band (Zong et al. 2018;Wang et al. 2019). Since the spectra were reduced, they were to be used directly. Before analysis, the spectra were normalized using Chebyshev function provided in the astropy package. An example of the normalized spectra is shown in Figure 12.
To get the RVs, the object spectra were matched with a template spectrum (spectral type K) to obtain CCF (cross correlation function) profiles and then applied double gaussian fit to the profiles. Only spectra in the blue band were used for CCF because the strong absorption of Hα in the red band may be easily affected by spot activities. The resulted CCF profiles are shown in Figure 13 and the fitting results are listed in Table 6, together with a brief log of observations (from the header of fits files). It should be noted that the errors in this table are only the errors of gaussian fit, so the real errors should be larger.
To compare the spectroscopic result with photometric solutions, the W-D program was utilized to calculate the RVs which takes consideration of eclipse-proximity correction (see the    Table 7 (the third row).
Here we also take a glance at the Hα line of our target from the spectra. Figure 15 shows how the broad Hα absorption changes with time. The absorption line from the primary component (more massive) seems to be missing or even has emission feature (msp119). This is peculiar and might indicate the strong activity on the surface of the binary.

Discussions and Conclusions
Comprehensive photometric solutions were derived independently both for ground-based data and space-based data. For the light curves from ground-based telescopes, a few cases including spots, a third light and even spot evolution were taken into consideration. For the light curves from K2 survey, cases with spots were mostly considered. The results from all cases indicate that the variable is a W-subtype shallow contact binary with a moderate high orbital inclination (i = 73.5 • ± 0.5 • ). However, the solutions vary in some parameters, especially the mass ratio q.
It is hard to figure out which is the correct one simply from the photometric solutions because even all of them give good fit (see Figure 7). It is suggested that the determination of mass ratio is affected by the spot settings. The spectroscopic investigation based on LAMOST msp data helps to clarify the result.

W-subtype Active Contact Binary System
Based on the results from both photometric solutions and RV analysis, the mass ratio of IL Cnc is calculated to be q = M 2 /M 1 = 1.76 ± 0.05(1/q = 0.57 ± 0.02). This value happens to be within the range that Alton (2018) declared. Together with the temperature difference about ∆T = T 2 − T 1 = −280 ± 20 K and a fill-out factor approximately 9%, it is clarified that IL Cnc is a W-subtype (the temperature of the less massive component is higher than that of the more massive component) shallow contact system. The RV curves give a further confirmation for this. This seems to support that K-type contact binary stars are more likely to be W-subtype.
The asymmetry in the light curves was usually modeled by spots on the components which indicates that the system may posses solar-like activities (e.g. Qian et al. (2007)). The continuous change of O'Connell effect in the light curves from K2 survey further tells that this system is highly active. The peculiar feature of Hα absorption of the spectral data may possibly be caused by spot activities. Another issue to be mentioned is that the errors listed in almost all the photometric solutions are pretty small, compared to the real uncertainties inferred from solutions of all cases. This suggests that errors from single set of solutions need to be adopted prudently.
Except for using the RV data as mentioned in Section 5, the absolute parameters of IL Cnc could also be estimated directly based on the mass-luminosity diagram, which is shown in Figure 16. Taking advantage of the parallax data from Gaia, the luminosity value was calculated with the following formulae (e.g. Chen et al. (2018)), where the parallax π = 4.08 ± 0.04 mas (Gaia Collaboration et al. 2018), m V = 12.6 mag (VSX database 6 ), the interstellar extinction A V ∼ 0 (b > 30 • ), and BC = −0.33 (Bessell et al. 1998).
The total luminosity L tot is then estimated to be about 0.579 L . By using the W-D code computed mass-luminosity line (based on the spot solutions from LCs 2016, see Table 4), two values were derived for the mass of IL Cnc based on whether it follows the classical massluminosity relation (hereafter CMLR, Eker et al. (2015)) or it exactly meets the estimation from Gaia parallax (see Figure 16). The estimated values of absolute parameters are listed in Table 7. Implied from Figure 16

Supplemental comments on the orbital period investigation
The results on orbital period investigation were presented both in Section 3.1 and 4.1, which are based on only groud-based data and the whole data, respectively. Besides the linear fit, we tried to fit the O-C diagrams with a general equation (Irwin 1952): 1+e cos ν sin(ν + w) + e sin w] Detailed explanation of the parameters can be found in the paper of Liao & Qian (2010). The red solid line (thick) denotes a classical mass-luminosity relation (CMLR). The green solid (thin) and dotted lines denote the total luminosity (Ltot) and the luminosity of star 2 (L2) that computed by W-D code, respectively. The blue dotted line represents the luminosity value based on Gaia parallax data. The cross denote the result from RV data.
A weighted least-squares method was used in the fitting. The results turned out that only the groud-based data could have a rough solution. After the K2 data were employed, no meaningful solution was found, which probably suggests that no eccentric LITE (light-time effect, see Frieboes-Conde & Herczeg (1973)) solution works here. Taking account that spot activities might affect the determination of the real conjunction of eclipse, we tried to fit the overall O-C diagram separately for the primary and secondary minima. Finally, the periodic solution was only found for the primary data. According to this solution, the orbital period changes with a period about 10.6 yr.
One explanation for the cyclic variation might be the light travel time effect (sometimes called "light-time effect"). Using a well known method (an example can be found in Liu et al. (2015a)) based on the periodic solution, the mass function was calculated to be f (m) = 8.8(±2.6) × 10 −5 M and the minimum mass of an additional body was then derived to be m 3 (i = 90 • ) = 0.061 ± 0.006 M . However, the difficulty of this hypothesis here is that the overall O-C diagram has a large scatter which could not be well fitted by the periodic equation and the cyclic variation only matches the primary minima.
Another hypothesis seems to be more convincing, which claims that the oscillation may be caused by the spot activity mechanism. Recently, Tran et al. (2013) reported that the O-C curves of short period binary stars display quasi-periodicities with typical amplitudes of 200∼300 s due to spot activities. Balaji et al. (2015) Figure 8 helps to demonstrate that the spot activity mechanism may be more plausible for the cyclic variation. It should be mentioned that on the opposite, Rappaport et al. (2013) reported a few candidate triple systems based on ETV (eclipsing time variation) data, a few of which have amplitudes below hundreds of seconds. This suggests that LITE explanation is not rejected yet, taking account of the relatively large errors (∼ 0.001 days) of the early timing data. To confirm the main reason of the cyclic variation here (tentatively spot activity mechanism), more precise data (especially continuous data) might be required in the future.
For secular variation of the orbital period, the results from the O−C diagrams suggest very small values (the period change rate (absolute value) is determined to be only 7(±5)×10 −9 days · yr −1 from Eq. 4). This is quite similar to some early K-type systems such as V1799 Ori and RV CVn (Liu et al. 2014a), which indicates that this phenomenon may be prevalent among K-type contact binaries.