apply more corrections

This commit is contained in:
2025-06-05 20:23:54 +02:00
parent 1135d65ab4
commit 4b0fcff9b5
11 changed files with 193 additions and 136 deletions
@@ -3,14 +3,14 @@
\chapter{Data and Methods}
\label{sec:data}
In this chapter the selection criteria for the data are explained in section \ref{sec:data:data_selection}. Furthermore a detailed description of the algorithm is given in section \ref{sec:data:data_reduction}.
In this chapter the selection criteria for the data are explained (section \ref{sec:data:data_selection}). Furthermore a detailed description of the algorithm is given in section \ref{sec:data:data_reduction}.
\section{Data selection}
\label{sec:data:data_selection}
This study uses a broad selection of Kepler/K2 and TESS lightcurves (427 for Kepler, 27 for K2, 1579 for TESS), which were downloaded from MAST with the help of the astroquery (\cite{astroquery}) python package. For Kepler and its continuation mission K2 short-cadence data were used. Unlike long-cadence data with a cadence of 30 minutes, this allows the resolution of shorter events too, as the duration of flares can vary between a few tens of seconds to a few hours (\cite{flare_duration1}).
The initial dataset was taken from a list of well known flaring stars from \cite{kepler_411_study}, \cite{kepler_411_210_comparison}, \cite{doyle_2018}, \cite{doyle_2019}, \cite{au_mic_flaring_spi} and \cite{flare_occurance_periodicity}. These stars were selected because previous studies showed that these stars show detectable flares. This allowed for a test of the algorithm described in the next section and compare it the previous results.
Additionally the dataset of the M to F stars of \cite{althukair_starlist} which can be found at \href{https://github.com/akthukair/AFD}{https://github.com/akthukair/AFD} was parsed to the downloader GUI. Due to not all data being available as short cadence data from Kepler/K2, only a subset of this large sample was added. The full list of stars, for which datasets were downloaded and parsed can be found in appendix \ref{chap:list_of_stars} table \ref{apA:list_of_all_stars}.
The initial dataset was taken from a list of well known flaring stars from \cite{kepler_411_study}, \cite{kepler_411_210_comparison}, \cite{doyle_2018}, \cite{doyle_2019}, \cite{au_mic_flaring_spi} and \cite{flare_occurance_periodicity}. These stars were selected because previous studies revealed that these stars showed detectable flares. This allowed a test of the algorithm described in the next section and a comparison with previous results.
Additionally the dataset of M to F stars of \cite{althukair_starlist} which can be found at \href{https://github.com/akthukair/AFD}{https://github.com/akthukair/AFD} was parsed to the downloader GUI. As not all data were available as short cadence data from Kepler/K2, only a subset of this large sample was added. The full list of stars, for which datasets were downloaded and parsed can be found in appendix \ref{chap:list_of_stars} table \ref{apA:list_of_all_stars}.
\section{Data reduction algorithms}
\label{sec:data:data_reduction}
@@ -24,32 +24,33 @@ The first step is to normalize the lightcurve. This is done to apply the same th
The next step is then to call $calculateFlareFitsForLightcurve()$ with the flattened lightcurve as well as the normalized lightcurve as parameters. It returns two lists of python dictionaries with the data for the flare peak as well as a fit which is described in the following paragraphs.
It parses the flattened lightcurve with the scipy $find\_peaks$ function. This function returns local maxima, which can be further filtered by their minimum height (optional) as well as the minimum distance of datapoints (optional). Afterwards the found peaks are sorted by height, and the highest 100 are returned. This was found to be a good amount as the most flares per fits file found in this study were around 70 for CD-56 1032A and B.\\
Afterwards each individual peak is checked. The first step is to estimate a flare starting and end point. This is done by checking the datapoints before and after the peak. If it finds that the flux delta is below 0.005 for three consecutive datapoints, it stops, and assumes that the last checked point is the start/end of the flare. In the case it finds an infinite or NaN value (which can happen if there are gaps in the lightcurve data), or it reaches 100 datapoints before/after it will stop. This was found to cover most flares detected and provides a sufficient number of datapoints for the following steps.\\
Afterwards multiple checks are done. It first checks if the datapoints right before and after the peak are above a threshold of 0.003 above the mean flattened flux, or if the two datapoints right after the peak are above that threshold (which is a similar approach to \cite{kepler_411_study}). Afterwards it is checked if the datapoint at two indices before the peak is larger than the datapoint right before the registered peak. This eliminates possible false positives which were found by manually inspecting previously found events, but could not be guaranteed to be two seperate events. In such cases it will only count one flare. Then a fit of the flare is generated. The first half of the fit, till the peak, is that of a gaussian function, with the second half being an exponential decay (similar approach to \cite{au_mic_flaring_spi} and \cite{doyle_2018}). Then the residual sum of squares (RSS) between the fit and the flux of the flare, as well as the total sum of squares (TSS) are calculated. Afterwards R-squared is calculated, and if it is below 0.8, the flare is rejected as the flare would not have the typical form. In the last step, the locations of the flare in the normalized and flattened lightcurve are compared. This step has been introduced, as in some rare cases the flattening algorithm can produce a large artificial spike (values of 10 or higher when normalized).
Afterwards multiple checks are done. It first checks if the datapoints right before and after the peak are above a threshold of 0.003 (0.3\%) above the mean flattened flux, or if the two datapoints right after the peak are above that threshold (which is a similar approach to \cite{kepler_411_study}). Afterwards it is checked if the datapoint at two indices before the peak is larger than the datapoint right before the registered peak. This eliminates possible false positives which were found by manually inspecting previously found events, but could not be guaranteed to be two seperate events. In such cases it will only count one flare. Then the flare lightcurve is fitted. For the impulsive phase, a gaussian function is used and for the gradual phase an exponetial function is used to account for the exponential decay in flare lightcurves (similar approach to \cite{au_mic_flaring_spi} and \cite{doyle_2018}). Then the residual sum of squares (RSS) between the fit and the flux of the flare, as well as the total sum of squares (TSS) are calculated. Afterwards R-squared is calculated, and if it is below 0.8, the flare is rejected as the flare would not have the typical form. In the last step, the locations of the flare in the normalized and flattened lightcurve are compared. This step has been introduced, as in some rare cases the flattening algorithm can produce a large artificial spike (values of 10 or higher when normalized).
\subsection{Lightcurve folding}
\label{sec:data:data_reduction:lightcurve_folding}
This section mainly describe how the $getOptimizedFold()$ function works. This function was written to optimize the lightcurve folding, to improve accuracy for the folding and lightcurve fitting process. It takes the normalized lightcurve as well as the preferred fitting function type (see chapter \ref{sec:gui:data_display}) as parameters. The fit type can either be "sine" for a sine fit, "poly" for a polynomlial fit, or "linear" for a linear fit. The default value is "sine", but can be changed for each individual star in the GUI.
The function at first generates two periodograms with the lightkurve function $to\_periodogram$. The first one uses the lombscargle algorithm, while the second one uses the boxleastsquares algorithm. Afterwards the 4 highest peaks of each are taken and converted into periods (unit in days).
This section mainly describes how the $getOptimizedFold()$ function works. This function was written to optimize the lightcurve folding, to improve accuracy for the folding and lightcurve fitting process. It uses the normalized lightcurve as well as the preferred fitting function type (see chapter \ref{sec:gui:data_display}) as parameters. The fit type can either be "sine" for a sine fit, "poly" for a polynomlial fit, or "linear" for a linear fit. The default value is "sine", but can be changed for each individual star in the GUI.
The function at first generates two periodograms with the lightkurve function $to\_periodogram$. The first one uses the lombscargle algorithm, while the second one uses the box-least-squares algorithm. Afterwards the four highest peaks of each are taken and converted into periods (unit in days).
Then each period found by each algorithm is compared with the periods found by the other, and in the case of a match (absolute value of the difference between the values of the two algorithms is smaller than 5\% of the larger of the two periods) this is now used as the rotational period as well as the spot modulation period. If it does not find a match, it uses the period corresponding to the highest peak in the lombscargle periodogram.
Afterwards it calculates the epoch time to the first minimum in the lightcurve using numpys $argrelextrema$ function. In the following loop, which is repeated up to 30 times, the lightcurve is folded (using the lightkurve $fold$ function) with the spot modulation period, and the epoch set to the first minima in the lightcurve. Then a fit is performed using the prefered fitting method. If the fitting process fails for the preferred fitting type, and it was set to "site", it will fall back to "poly" (polynomial fit), and vice versa. If the resulting fit has two maxima, and they are further away from the edges than 10\% of the used period, the algorithm checks if there is a signal for half the used period in any of the two periodograms. If this is the case, it will use this as the period for spot modulation. It will save the current folded lightcurve as folded by rotational period, and generate a new folded lightcurve with the newly found spot modulation period. Afterwards it checks the location of the minimum of the spot modulation folded lightcurve. If the minimum is within 1\% of half the phase it will stop. If not, it will shift the epoch by the difference of the minimum to 0, which is the center of the folded lightcurves phase, and repeat the folding process.\\
The function then returns the rotational period, the spot modulation period, the periods found by both periodograms as well as the folded lightcurve, the corresponding phase and the fit and fit type ("sine", "poly" or "linear") as well as the epoch. If the spot modulation differs from the rotational period found, it will also return the folded lightcurve, phase, fit and fit type of the rotational period folded lightcurve. The last returned value is "isValid", which is set to false if the algorithm does not find a fit with its minimum at half the phase for the spot modulation period folded lightcurve within 30 iterations.
Afterwards it calculates the epoch time to the first minimum in the lightcurve using numpys $argrelextrema$ function. In the following loop, which is repeated up to 30 times, the lightcurve is folded (using the lightkurve $fold$ function) with the spot modulation period, and the epoch set to the first minima in the lightcurve. Then a fit is performed using the prefered fitting method. If the fitting process fails for the preferred fitting type, and it was set to "sine", it will fall back to "poly" (polynomial fit), and vice versa. If the resulting fit has two maxima, and they are further away from the edges than 10\% of the used period, the algorithm checks if there is a signal for half the used period in any of the two periodograms. If this is the case, it will use this as the period for spot modulation. It will save the current folded lightcurve as folded by rotational period, and generate a new folded lightcurve with the newly found spot modulation period. Afterwards it checks the location of the minimum of the spot modulation folded lightcurve. If the minimum is within 1\% of half the phase it will stop. If not, it will shift the epoch by the difference of the minimum to zero, which is the center of the folded lightcurve phase, and repeat the folding process.\\
The function then returns the rotational period, the spot modulation period, the periods found by both periodograms as well as the folded lightcurve, the corresponding phase and the fit and fit type ("sine", "poly" or "linear") as well as the epoch. If the spot modulation differs from the found rotational period, it will also return the folded lightcurve, phase, fit and fit type of the rotational period folded lightcurve. The last returned value is "isValid", which is set to false if the algorithm does not find a fit with its minimum at half the phase for the spot modulation period folded lightcurve within 30 iterations.
\section{Spectral type identification}
\label{sec:data:sptype_identification}
After the initial download and cross checking with SIMBAD, not all stars had a spectral type assigned. For these stars magnitudes in different wavelength bands was available. The stars with known magnitudes known atleast in the B and V bands are given in table \ref{tab:unknown_sptypes_bv}, while the stars with only J, H and K bands are in table \ref{tab:unknown_sptypes_jhk}.
After the initial download and cross checking with SIMBAD, not all stars had a spectral type assigned. For these stars magnitudes in different wavelength bands were available. The stars with known magnitudes at least in the B and V bands are given in table \ref{tab:unknown_sptypes_bv}, while the stars with only J, H and K bands are given in table \ref{tab:unknown_sptypes_jhk}.
Using $B-V$ (as well as further magnitude differences, if available) have been compared to the table at \cite{spectral_type_color_table} (\href{https://www.stsci.edu/~inr/intrins.html}{https://www.stsci.edu/~inr/intrins.html}).
The B-V color index as well as other color indices (if available) have been compared to the table in \cite{spectral_type_color_table} (\href{https://www.stsci.edu/~inr/intrins.html}{https://www.stsci.edu/~inr/intrins.html}).
This allowed a rough classification of the stars in table \ref{tab:unknown_sptypes_bv}. Additionally, also the effective temperature was taken into account from the TESS Input Catalogue (\cite{revised_tess_input_catalogue}) from Vizier (\cite{vizier}). The classification table for spectral types based on effective temperature has been taken from \cite{harvard_spectral_types_teff} (\href{https://lweb.cfa.harvard.edu/~pberlind/atlas/htmls/note.html}{https://lweb.cfa.harvard.edu/~pberlind/atlas/htmls/note.html}).\\
The final spectral type for the stars can be found in table \ref{tab:unknown_sptypes_final_output}.
If the color spectral type calculated from the effective temperature and the magnitude differences in multiple different color bands is identical, it will use that spectral type as final.
If the spectral type calculated from the effective temperature and the indices is identical, it will use that spectral type as final.
In the case of those differing, like for BD-08 995, which was found to be an early K type star (K0) by the intrinsic color method but being a G type star based on the effective temperature, the latter is prefered. For the four stars from table \ref{tab:unknown_sptypes_jhk} which have no B or V measurements, the spectral type based on the effective temperature is taken as final.
In the case of those differing, like for BD-08 995, which was found to be an early K type star (K0) by the intrinsic color method but being a G type star based on the effective temperature, the latter is prefered. For the four stars from table \ref{tab:unknown_sptypes_jhk} which have no B or V measurements, the spectral type based on the effective temperature is used.
\begin{table}
\centering
\caption{List of stars with no spectral type entry in SIMBAD, but available magnitudes in the B and V bands. Main Identifier according to SIMBADs "MAIN\_ID" property, as well as TESS Input Catalogue and Kepler Input Catalogue numbers are given.}
\label{tab:unknown_sptypes_bv}
\begin{tabular}{lll}
@@ -68,7 +69,8 @@ In the case of those differing, like for BD-08 995, which was found to be an ear
\begin{table}
\caption{List of stars with no spectral type entry on SIMBAD, but available magnitudes in the J, H and K bands. Main Identifier according to SIMBADs "MAIN\_ID" property (2MASS), as well as TESS Input Catalogue and Kepler Input Catalogue numbers are given.}
\centering
\caption{List of stars with no spectral type entry in SIMBAD, but available magnitudes in the J, H and K bands. Main Identifier according to SIMBADs "MAIN\_ID" property (2MASS), as well as TESS Input Catalogue and Kepler Input Catalogue numbers are given.}
\label{tab:unknown_sptypes_jhk}
\begin{tabular}{lll}
\hline
@@ -84,6 +86,7 @@ In the case of those differing, like for BD-08 995, which was found to be an ear
\begin{landscape}
\begin{table}
\centering
\caption{Final determined spectral types. Main Identifier taken from SIMBAD. The column B-V contains the value of B-V calculated from the respective magnitudes taken from SIMBAD. The columns T\textsubscript{eff}, MASS and Radius contain the effective temperature, mass and radius for the star taken from the TESS Input Catalogue from Vizier. The column Spectral Type (B-V) contains the determined spectral type based on \cite{spectral_type_color_table}, while the column Spectral Type (T\textsubscript{eff}) contains the spectral type based on the effective temperature from \cite{harvard_spectral_types_teff}. The last column contains the spectral type which is used in this study.}
\label{tab:unknown_sptypes_final_output}
\begin{tabular}{lccccccc}