The VARTOOLS Light Curve Analysis Program


The VARTOOLS program is a command line utility that provides tools for calculating variability/periodicity statistics of light curves as well as tools for modifying light curves. The program is run by issuing a sequence of commands to perform actions on light curves, each command is executed in turn with the resulting light curves passed to the next command. Statistics computed by each command are sent to stdout as an ascii table.

This program is available without any warranty - use it at your own risk. If you do find a bug, please let me know so that I can fix it (jhartman AT cfa DOT harvard DOT edu). Also feel free to let me know if there are any commands/algorithms that you would like to see added to this program.

Citation: If you use this program, please cite Hartman et al., 2008, ApJ, 675, 1254. Also, be sure to cite the relavent sources for the algorithms that you use.

Download

To install, unpack the program to its own directory. Then type "make" from that directory. Type "./vartools" for the syntax, or "./vartools -help" for a detailed description. The program is written in ansi c and only uses standard libraries, so it should compile on any linux/unix system without difficulty.

Commands: In the current version, the following commands are supported:

-alarm :

Calculate the alarm variability statistic for each light curve. Cite Tamuz, Mazeh, and North 2006, MNRAS, 367, 1521 if you use this tool.

-aov ["Nbin" Nbin] minp maxp subsample finetune Npeaks operiodogram [outdir] :

Perform an AoV period search on the light curves. This implementation uses phase-binning. Specify "Nbin" and a number to change the number of bins used from the default value of 8. minp and maxp are the minimum and maximum periods to search. The intial search will use a frequency resolution of subsample/T where T is the baseline of the light cruve. The peak periods will be refined using a resolution of finetune/T. The program will output the Npeaks highest peaks in the periodogram of the light curve. The output aov statistic for each light curve is (<-ln(theta_aov)> - ln(theta_aov))/RMS(-ln(theta_aov)) where the average and RMS are taken over the entire periodogram for a light curve, the output also includes the average and RMS of -ln(theta_aov). operiodogram should be either 0 or 1. If it is set to 1 then the periodogram for each light curve will also be output to the directory outdir, with the suffix ".aov". The first column in the output is the period, the second column is theta_aov. Cite Schwarzenberg-Czerny, A. 1989, MNRAS, 241, 153 and Devor, J. 2005, ApJ, 628, 411 if you use this tool. (The Devor citation is needed because this code uses his implementation of AoV).

-aov_harm Nharm minp maxp subsample finetune Npeaks operiodogram [outdir] :

Perform an AoV period search on the light curves using a multi-harmonic model. Nharm is the number of harmonics to use. minp and maxp are the minimum and maximum periods to search. The initial search will use a frequency resolution of subsample/T where T is the baseline of the light cruve. The peak periods will be refined using a resolution of finetune/T. The program will output the Npeaks highest peaks in the periodogram of the light curve. The output statistic for each light curve is theta_aov (note this is different from the normalized -ln(theta_aov) that is returned by the -aov command), the output also includes the average and RMS of theta_aov. operiodogram should be either 0 or 1. If it is set to 1 then the periodogram for each light curve will also be output to the directory outdir, with the suffix ".aov_harm". The first column in the output is the period, the second column is theta_aov. Cite Schwarzenberg-Czerny, A., 1996, ApJ, 460, L107 if you use this tool.

-binlc combinetype <"binsize" binsize | "nbins" nbins> ["firstbinshift" firstbinshift] timetype :

Bin the light curves in time (or phase if a -Phase command has already been issued). combinetype is 0 to take the average of points in a bin, 1 to take the median of the points or 2 to take the weighted average. One should specify either the binsize in units of the time coordinate (e.g. in phase if the light curves have been phased), or the number of bins to split the light curve time-span into. By default the first bin begins at the initial time in the light curve (t0), if firstbinshift is specified, then this will be shifted by t0 - firstbinshift/binsize. timetype should be 0, 1 or 2. If it is 0 then the output time for each bin is the time at the center of the bin, if it is 1 then the time will be the average of the times of points that fall within the bin, if it is 2 then it will be the median.

-BLS < "r" rmin rmax | "q" qmin qmax > minper maxper nfreq nbins timezone Npeak outperiodogram [outdir] omodel [modeloutdir] correctlc :

This command runs BLS on the light curves. You may either set a fixed minimum and maximum q (fraction of orbit in transit) for the search, or you can specify a minimum and maximum stellar radius to consider (in solar radii) in which case the qmin and qmax values are calculated for each trial period P as q = 0.076 * R**(2/3) / P**(2/3) (with P in days - this assumes R = M as an approximation for the lower main sequence). The minimum and maximum periods to search are given in days, nfreq is the number of trial frequencies (10000 is a decent number). nbins is the number of phase bins to break the light curve into (200 is a decent number). timezone is the number to add to the julian date to get the local date (-7 for arizona), this is used to determine which nights different observations come from and thus to determine the fraction of delta-chi2 that comes from a single night. Npeak is the number of peaks in the BLS spectrum to find and report. outperiodogram is a flag that is set to 1 to output the BLS period vs. SN spectrum. outdir is the output directory for the BLS spectrum if outperiodogram is set to 1, ".bls" will be appended to the filename. omodel is a flag set to 1 or zero that can be used to output the model for the light curve, the output directory is then given in modeloutdir, the suffix ".bls.model" will be appended to the filename. correctlc is either 0 or 1, set to 1 it will subtract the transit model from the light curve before passing it to the next command. Cite Kovács, G., Zucker, S., & Mazeh, T. 2002, A&A, 391, 369 if you use this tool.

-BLSFixPer <"aov" | "ls" | "spec"> <"r" rmin rmax | "q" qmin qmax > nbins timezone omodel [model_outdir] correctlc :

This command runs BLS at a fixed period on the light curves - it searches merely for the most transit-like signal at the specified period. The period comes either from the last aov command, from the last ls command, or is specified (in which case it must be in the appropriate column of the lc_list). Like BLS you can either specify a minimum and maximum stellar radius to consider, or a minimum and maximum q to consider. Similarly you must specify the number of bins to use and the timezone. omodel is 1 to output the model light curve to model_outdir (with ".blsfixper.model" appended to the end of the file name. correctlc is 1 to subtract the transit model from the light curve before passing it to the next command. Cite Kovács, G., Zucker, S., & Mazeh, T. 2002, A&A, 391, 369 if you use this tool.

-chi2 :

Calculate chi2 per dof for the light curves. The output will include chi2 and the error weighted mean magnitude.

-chi2bin Nbin bintime1...bintimeN :

Calculate chi2 per dof after applying a moving mean filter to the light curves. (Note that the light curves passed to the next command are unchanged). Nbin filters are used (with Nbin resulting estimates of chi2 and the error weighted mean). The width of each filter is given by 2.0*bintime (bintime is given in minutes)

-clip sigclip iter :

Sigma-clip the light curves using a clipping factor of sigclip. iter is a flag that is either 1 for iterative clipping, or 0 for single pass clipping. The output table will include the number of points that were clipped.

-decorr correctlc zeropointterm subtractfirstterm Nglobalterms globalfile1 order1 ... Nlcterms lccolumn1 lcorder1 ... omodel [modeloutdir] :

Decorrelate the light curves against specified signals. The signals are either global signals (e.g. the airmass) that are input in separate files (with the format: JD signal_value), or are light curve specific signals (e.g. the sub-pixel position) that are input as additional columns in the light curves. correctlc is a flag that is 0 or 1, if the value is 1 then the light curves passed to the next command will be decorrelated, if the value is 0 then resulting chi2 from decorrelation and the decorrelation coefficients will be output to the table, but the light curves themselves will not be decorrelated. zeropointterm is a flag that is 0 or 1, if the value is 1 then a zeropoint term is included in the decorrelation, if it is 0 then no such term is included. subtractfirstterm is a flag that is 0 or 1, if the value is 1 then the light curves are decorrelated against the signal minus the first signal value (this is useful if one is decorrelating against the JD, for example, in which case one would decorrelate against JD - JD0 where JD0 is the first JD in the light curve), if it is 0 then the light curves are decorrelated against the signal without subtracting the first term. Nglobalterms is the number of global files. globalfile1...globalfileN or the names of those files, order1...orderN are the orders of the polynomials used to decorrelate each signal (each must be greater than 0). Nlcterms is the number of light curve specific signals. The columns of these signals are given by lccolumn1...lccolumnN. The orders of the polynomials are given by lcorder1...lcorderN. omodel is a flag set to 1 or zero that can be used to output the model for the light curve, the output directory is then given in modeloutdir, the suffix ".decorr.model" will be appended to the filename.

-ensemblerescalesig sigclip :

Rescale the magnitude uncertainties by fitting a linear relation between (expected rms)^2 vs (chi2)*(expected rms)^2 for all light curves so that chi2 per dof is distributed about unity, Outliers in the distribution are clipped using a clipping factor of sigclip. This routine requires a list of light curves. If this command is invoked, then all light curves will be read into memory. The output will include the average rescale factor for each light curve (taken to be sqrt(chi2_after/chi2_before)).

-fluxtomag mag_constant offset :

Convert light curves from ISIS differential flux into magnitudes. In this case a light curve list must be used, and the reference magnitudes of the stars (after aperture correction) should be given as an additional column in the light curve list. mag_constant is the magnitude of a source with a flux of 1 ADU, offset is an additive constant to apply to the output light curves. This command does not yield any output to stdout, you must explicitly call -rms or -chi2 if you want statistics.

-GetLSAmpThresh <"ls" | "spec"> minp thresh Nharm Nsubharm :

This routine determines the minimum amplitude that a light curve could have and still be detected at a given period (which is either specified in the input list, or taken from the last LS command) with a LS log(formal false alarm probability) < thresh. This is done by fitting a fourier series with Nharm harmonics and Nsubharm subharmonics to the light curve, subtracting it from the light curve and re-adding it after scaling it by a number. minP is the minimum period that would be searched (this is needed since it sets the false alarm probability scale). The output is the minimum factor by which the signal could be scaled and still be detectable together with the minimum peak-to-peak amplitude.

-Jstet timescale dates :

Calculate Stetson's J statistic, L statistic and the Kurtosis of each light curve. The timescale is the time in minutes that distinguishes between "near" and "far" observations. The dates file should contain JDs for all possible observations in the first columns - this is used to calculate the maximum possible weight. Cite Stetson, P.B. 1996, PASP, 108, 851 if you use this tool.

-Killharm <"aov" | "ls" | "both" | "spec" Nper> Nharm Nsubharm omodel [modeloutdir] :

This command whitens light curves against one or more periods. The mean value of the light curve and the cos and sin coefficients are output. The light curves passed to the next command are whitened light curves. The origin of the period(s) is either from the most recent previous aov command, from the most recent previous LS command, two periods one from aov, the other from LS, or Nper periods specified in the input list (the periods are read off in order as additional columns in the input light curve list - a list must be used for this option). The light curves will be whiten using Nharm higher-harmonics and Nsubharm sub-harmonics. omodel is a flag set to 1 or zero that can be used to output the model for the light curve, the output directory is then given in modeloutdir, the suffix ".killharm.model" will be appended to the filename.

-LS minp maxp subsample Npeaks operiodogram [outdir] :

Perform a Lomb-Scargle period search on the light curves. In this case the statistic that is output is the false alarm probability for the peaks (assuming white noise with rms equal to the rms of the light curve). The output periodogram has the format: frequency normalized-LS statistic. This version does not allow the mean to float. Cite Lomb, N.R. 1976, A&SS, 39, 447, Scargle, J.D. 1982, ApJ, 263, 835, Press, W.H. & Rybicki, G.B. 1989, ApJ, 338, 277, and Press, W.H., Teukolsky, S.A., Vetterling, W.T. & Flannery, B.P. 1992, Numerical Recipes in C, 2nd ed. (New York: Cambridge University Press) if you use this tool.

-MandelAgolTransit <bls | blsfixper | P0 T00 r0 a0 sin_i0 e0 omega0 mconst0> <"quad" | "nonlin"> ldcoeff1_0 ... ldcoeffn_0 fitP fitT0 fitr fita fitsin_i fite fitomega fitmconst fitldcoeff1 ... fitldcoeffn correctlc omodel [model_outdir] :

Fit a Mandel and Agol (2002) transit model to the light curve. The initial parameters either are automatically estimated using the results from bls if "bls" is specified, a blsfixper command, or they are entered directly in the command line as the values for P0 T00 etc. The parameters are P0 - period, T00 - time of transit, r0 - ratio of planet radius to star radius, a0 - ratio of semi-major axis to star radius, sin_i0 - sine of the inclination angle, e0 - eccentricity, omega0 - argument of periastron, and mconst0 - the out of transit magnitude. If mconst0 is negative then it will be estimated directly from the light curve. "quad" or "nonlin" are used to specify the limb darkening law, either quadratic or non-linear, following Claret. The number of initial limb darkening coefficients must be 2 for quadratic limb-darkening or 4 for non-linear limb-darkening. fitP fitT0 ... should be 1 or 0 depending on whether or not the corresponding parameter is allowed to vary. Set correctlc to 1 to subtract the model from the light curve, set it to 0 to only perform the fit. omodel is a flag set to 1 or zero that can be used to output the model for the light curve, the output directory is then given in modeloutdir, the suffix ".mandelagoltransit.model" will be appended to the filename. Cite Mandel, K., & Agol, E. 2002, ApJ, 580, L171 if you use this tool.

-o outdir :

output the light curves to directory outdir

-Phase <aov | ls | bls | spec> :

Replace the time variable of a light curve with its phase and sort the light curve by the phase. The period for phasing is either taken from the last aov command, the last ls command, the last BLS command, or is specified in the list (in this case the list must be used).

-rescalesig :

Rescale the magnitude uncertainties of each light curve so that Chi2 per dof is equal to 1 for every light curve. The rescale factor for each light curve will be included in the output table.

-rms :

Calculate the rms of the light curves. The output will include RMS, the mean magnitude, the expected RMS and the number of points in the light curve

-rmsbin Nbin bintime1...bintimeN :

Similar to chi2bin, this calculates RMS after applying a moving mean filter.

-SoftenedTransit <bls | blsfixper | P0 T00 eta0 delta0 mconst0> cval0 fitP fitT0 fiteta fitcval fitdelta fitmconst correctlc omodel [model_outdir] fit_harm [<"aov" | "ls" | "bls" | "spec" Pharm> nharm nsubharm] :

Fit a Protopapas, Jimenez and Alcock transit model to the light curve. The initial parameters either come from bls if "bls" is specified, a blsfixper command if "blsfixper" is specified, otherwise they are entered directly in the command line as the values for P0 T0 etc. If mconst0 is < 0 then it will be estimated directly from the light curve. fitP fitT0 should be 1 or 0 depending on whether or not the corresponding parameter is allowed to vary. Set correctlc to 1 to subtract the model from the light curve, set it to 0 to only perform the fit. omodel is a flag set to 1 or zero that can be used to output the model for the light curve, the output directory is then given in modeloutdir, the suffix ".softenedtransit.model" will be appended to the filename. If fit_harm is 1 then simultaneously fit a harmonic series to the light curve. The period comes from either the last aov, the last ls, the last bls model or is specified in the command line as Pharm (you must type "spec" and then the value), finally the number of harmonics and sub-harmonics must also be specified. Cite Protopapas, P., Jimenez, R., & Alcock, C. 2005, MNRAS, 362, 460 if you use this tool.

-Starspot <aov | ls | spec> a0 b0 alpha0 i0 chi0 psi00 mconst0 fitP fita fitb fitalpha fiti fitchi fitpsi0 fitmconst correctlc omodel [modeloutdir] :

This command fits a single, circular, uniform temperature spot model to the light curve using the Dorren 1987 model. Use aov or ls to choose the initial period to be from AoV or LS, use spec to specify the initial period in the input list (the period is read off as an additional column in the input light curve list - a list must be used for this option). The other parameters are initial guesses as in Dorren 1987, mconst is the constant magnitude term. Set mconst0 to a negative value to have the program automatically estimate the initial mconst. Note that here we take a = a_d/(pi * (1 - mu_* / 3)) and b = b_d/(pi * (1 - mu_* / 3)) where a_d and b_d are the a and b terms from Dorren 1987. fitP... are flags denoting whether or not a parameter is to be varied. Set correctlc to 1 to subtract the model from the light curve, set it to 0 to only perform the fit. omodel is a flag set to 1 or zero that can be used to output the model for the light curve, the output directory is then given in modeloutdir, the suffix ".starspot.model" will be appended to the filename. Cite Dorren 1987, ApJ, 320, 756 if you use this tool.

-SYSREM Ninput_color Ninput_airmass initial_airmass_file sigma_clip1 sigma_clip2 saturation correctlc omodel [model_outdir] otrends [trend_outfile] useweights :

Run the SYSREM algorithm to identify and remove trends from an ensemble of light curves. If this command is used, the -readall option will automatically be set. Also, if this command is used then one must have as input a light curve list. A total of Ninput_color + Ninput_airmass trends will be searched for. Ninput_color of these will have the colors set initially, Ninput_airmass will have the airmass terms set initially. The initial color terms will be read in as additional columns in the input light curve list, there must be Ninput_color of these columns. The initial_airmass_file gives a list of the initial airmass trends to use, the first column is JD and the subsequent Ninput_airmass columns are the initial airmass trends. As implemented this routine first filters the trends with the airmass terms specified initially and then filters the trends with the color terms specified initially. sigma_clip1 is the sigma clipping used in calculating the mean magnitudes for the light curves. sigma_clip2 is the sigma clipping used in determining whether or not points contribute to the airmass or color terms when doing the fit. Any points with magnitude less than saturation will not contribute. omodel is a flag that is 1 or zero for outputing the model light curves, the output directory is given in model_outdir, the suffix ".sysrem.model" will be appended to the filename. The format for the model light curve is: JD mag mag_model sig clip, where clip is either 1 if the point was included in determining the trends or 0 if it wasn't. otrends is a flag that is 1 or zero for outputing the final trends. These will be output to the file trend_outfile, the first column is the JD and subsequent columns are for each trend signal. If you use this command, be sure to cite Tamuz, Mazeh and Zucker (2005, MNRAS, 356, 1466).

-TFA trendlist dates_file pixelsep correctlc ocoeff [coeff_outdir] omodel [model_outdir] :

Run the Trend Filtering Algorithm on the light curves. The trendlist is a file containing a list of light curves to be used as trends, it has the format: trendname trendx trendy, where trendname is the file name and trendx and trendy are the coordinates of the trend star. dates_file is a file giving all the dates in the second column, the first column is ignored (it is the same format as the dates file for the ISIS image subtraction program where the first column is the filename and the second is JD). Trend stars within pixelsep of the light curve in question will not be used in detrending the light curve. To you this routine you need an input light curve list, the x and y positions of each light curve must be given as additional columns in the list. correctlc is a flag that is 0 or 1 indicating whether the light curves passed to the next command should have the tfa applied. ocoeff is a flag to denote whether or not to output the list of coefficients for each trend for a given light curve, if set to 1 then they will be output to coeff_outdir with ".tfa.coeff" appended to the end of the filename. omodel is a flag set to 1 or zero that can be used to output the tfa model for the light curve, the output directory is then given in model_outdir, the suffix ".tfa.model" will be appended to the filename. If you use this routine be sure to cite Kovacs, Bakos and Noyes, 2005, MNRAS, 356, 557

Options: In the current version, the following options are available for specifying the input/output format:

<-i lcname | -l lclist> :

either provide an individual light curve, or a list of light curves. If lcname or lclist is "-" then input will be taken from stdin. Some commands can only be used with lists. The list should contain a single light curve filename per line in the first column. Additional columns may be necessary to specify, for example, the periods to use for pre-whitening or any other values required by some commands.

-header :

use this option to provide a one line header at the start of the output table.

-tab :

use this option to use a tab-delimited starbase format for the output table rather than the default space-delimited format.

-basename :

use this option to print only the basename of each light curve in the output table.

-readformat Nskip col_jd col_mag col_sig :

use this option to specify the format of the input light curves. Nskip is the number of lines to skip at the beginning of each file (this includes any lines that begin with '#' which are otherwise automatically ignored), the default value is 0. col_jd, col_mag and col_sig are the columns that contain the time, magnitude and magnitude uncertainties, the default values are 1, 2 and 3. (Note that the time is assumed to be in days).

-readall :

use this option to force the program to read in all the light curves at once. If not specified, this mode will only be used if a command that requires storing all the light curves in memory is issued.

-ascii :

If this is specified then periodograms will be printed out in ascii format. Otherwise they will be stored in binary.

-jdtol jdtol :

Time measurements within jdtol of each other are considered equal. The default value is 0.000010 days.

Examples: The following are examples for commonly called tasks, the examples may be called from the VARTOOLS installation directory, all files referenced in these examples are included in the .tar.gz package and unpack to the EXAMPLES subdirectory. Lines after the '$' are commands that are typed at the prompt.

Example 1. Fitting a quadratic polynomial in JD to a list of light curves.
$ ./vartools -l EXAMPLES/lc_list -rms -decorr 1 1 1 0 1 1 2 0 -rms -tab

Name Mean_Mag_0 RMS_0 Expected_RMS_0 Npoints_0 Decorr_constant_term_1 Decorr_constant_term_err_1 \
LCColumn_1_coeff_1_1 LCColumn_1_coeff_err_1_1 LCColumn_1_coeff_2_1 LCColumn_1_coeff_err_2_1 Decorr_chi2_1 \
Mean_Mag_2 RMS_2 Expected_RMS_2  Npoints_2
---- ---------- ----- -------------- --------- ---------------------- -------------------------- \
-------------------- ------------------------ -------------------- ------------------------ ------------- \
---------- ----- --------------  ---------
EXAMPLES/1 10.24745 0.15944 0.00101 3122 10.0830375984825 0.0000325849746 0.0097933162509 0.0000059117875 \
0.0002554062775 0.0000001956124 6.68601 10.24728 0.00211 0.00101 3122
...

Fit a quadratic polynomial to the light curves given in the file EXAMPLES/lc_list. To do this we use no global terms and 1 lc term. For the lc term we use the first column in the light curve (the JD) and fit a second order polynomial in this term to each light curve. We fit for the zero-point term and correct the light curve (so that commands after -decorr will receive light curves with the best-fit quadratic polynomial removed, note that when the light curve is corrected the mean is kept constant), we also subtract the first term in the signal that we are decorrelating against (use JD - JD0 rather than JD, since JD*JD runs into round-off problems whereas (JD - JD0)*(JD - JD0) does not), but we do not output the corrected light curves to the disk. The rms is determined before and after the fit, the -tab option causes the output table to be in starbase format with headers for each column. To interpret the output, note that for light curve 1 we find that the best-fit quadratic has the form: 0.0002554062775*(JD - 53725.173920)*(JD - 53725.173920) + 0.0097933162509*(JD - 53725.173920) + 10.0830375984825, and that fitting this equation reduces the RMS from 0.15944 mag to 0.00211 mag (a quadratic signal was injected into this particular light curve).

Example 2. Performing a Lomb-Scargle period search on a light curve and fitting a harmonic series to the light curve.
$ ./vartools -i EXAMPLES/2 -LS 1.0 2.0 0.01 1 0 -Killharm ls 0 0 1 EXAMPLES/OUTDIR1 -tab

Name LS_Period_1_0 Log10_LS_Prob_1_0 Mean_Mag_1 LS_Fundamental_Sincoeff_1 LS_Fundamental_Coscoeff_1 \
LS_Amplitude_1
---- ------------- ----------------- ---------- ------------------------- ------------------------- \
--------------
EXAMPLES/2 1.23588006 -707.23302 10.11178 0.02004 -0.04580 0.09998

Run the Lomb-Scargle period search on the light curve EXAMPLES/2. This light curve has had a sine-curve with a period of 1.2354 days and an amplitude of 0.05 mag injected into it. After using L-S to find the period, fit a harmonic series with no harmonics or sub-harmonics (i.e., just a sine-curve at the period found by L-S) to the light curve. Output the light curve with the best-fit model in the third column to the file EXAMPLES/OUTDIR1/2.killharm.model. Note that in this case L-S finds the period of 1.23588006 days, and the killharm routine fits the function 0.02004*sin((JD - JD0)*2*pi/1.23588006) - 0.0458*cos((JD - JD0)*2*pi/1.23588006) + 10.11178 to the light curve. The killharm routine returns a peak-to-peak amplitude of 0.09998 mag.

Example 3. Injecting a transit signal into a light curve.
$ ./vartools -i EXAMPLES/3 -MandelAgolTransit 2.12345 53725.174 0.1 10.0 1.0 0. 0. 0 quad 0.236 0.391 \
  0 0 0 0 0 0 0 0 0 0 0 1 EXAMPLES/OUTDIR1 -tab

Name MandelAgolTransit_Period_0 MandelAgolTransit_T0_0 MandelAgolTransit_r_0 MandelAgolTransit_a_0 \
MandelAgolTransit_sin_i_0 MandelAgolTransit_e_0 MandelAgolTransit_omega_0 MandelAgolTransit_mconst_0 \
MandelAgolTransit_ldcoeff1_0 MandelAgolTransit_ldcoeff2_0 MandelAgolTransit_chi2_0
---- -------------------------- ---------------------- --------------------- --------------------- \
------------------------- --------------------- ------------------------- -------------------------- \
---------------------------- ---------------------------- ------------------------
EXAMPLES/3 2.12345000 53725.17400000 0.10000 10.00000 1.00000 0.00000 0.00000 0.00000 0.23600 \
0.39100 129730808.07162

$ gawk '{print $1, $2 + $3, $4}' EXAMPLES/OUTDIR1/3.mandelagoltransit.model > EXAMPLES/3.transit

These two commands effectively inject a Mandel and Agol transit into a light curve. The first line calls the MandelAgolTransit vartools command to create a simulated transit with a period of 2.12345 days, T0 = 53725.174, Rplanet/Rstar = 0.1, a/Rstar = 10.0, sini = 1.0, a circular orbit, and quadratic limb darkening for the star with coefficients of 0.236 and 0.391. We do not fit for any of the parameters (10 zeros), nor do we correct the light curve (the 11th 0). The model is output in the third column of EXAMPLES/OUTDIR1/3.mandelagoltransit.model. In the second line we use gawk (or awk) to add the model transit signal to the light curve and output the result to the file EXAMPLES/3.transit.

Example 4. De-trending a light curve with TFA and running BLS on the de-trended light curve.
$ ./vartools -l EXAMPLES/lc_list_tfa -rms -TFA EXAMPLES/trendlist_tfa EXAMPLES/dates_tfa 25.0 1 0 0 \
  -BLS q 0.01 0.1 0.5 5.0 2000 200 7 2 0 1 EXAMPLES/OUTDIR1 1 -rms -tab

Name Mean_Mag_0 RMS_0 Expected_RMS_0 Npoints_0 TFA_MeanMag_1 TFA_RMS_1 BLS_Period_1_2 BLS_SN_1_2 BLS_SR_1_2 \
BLS_SDE_1_2 BLS_Depth_1_2 BLS_Qtran_1_2 BLS_i1_1_2 BLS_i2_1_2 BLS_deltaChi2_1_2 BLS_fraconenight_1_2 \
BLS_Npointsintransit_1_2 BLS_Ntransits_1_2 BLS_Npointsbeforetransit_1_2 BLS_Npointsaftertransit_1_2 \
BLS_Rednoise_1_2 BLS_Whitenoise_1_2 BLS_SignaltoPinknoise_1_2 BLS_Period_2_2 BLS_SN_2_2 BLS_SR_2_2 \
BLS_SDE_2_2 BLS_Depth_2_2 BLS_Qtran_2_2 BLS_i1_2_2 BLS_i2_2_2 BLS_deltaChi2_2_2 BLS_fraconenight_2_2 \
BLS_Npointsintransit_2_2 BLS_Ntransits_2_2 BLS_Npointsbeforetransit_2_2 BLS_Npointsaftertransit_2_2 \
BLS_Rednoise_2_2 BLS_Whitenoise_2_2 BLS_SignaltoPinknoise_2_2 BLS_Period_invtransit_3_2 \
BLS_deltaChi2_invtransit_3_2 BLS_MeanMag_3_2 Mean_Mag_3 RMS_3 Expected_RMS_3 Npoints_3
----  ---------- ----- -------------- --------- ------------- --------- -------------- ---------- ---------- \
----------- ------------- ------------- ---------- ---------- ----------------- -------------------- \
----------------------- ----------------- --------------------------- --------------------------- \
---------------- ------------------ ------------------------- -------------- ---------- ---------- \
----------- ------------- ------------- ---------- ---------- ----------------- -------------------- \
----------------------- ----------------- --------------------------- --------------------------- \
---------------- ------------------ ------------------------- ------------------------- \
---------------------------- --------------- ---------- ----- -------------- ---------
EXAMPLES/3.transit 10.16727 0.00542 0.00104 3417 10.16714 0.00471 2.12298216 72.07052 0.00243 4.84830 \
0.01162 0.03000 0.98500 0.01000 -25264.50212 0.42585 146 4 93 87 0.00140 0.00407 14.99298 1.35250338 \
59.87855 0.00204 3.55079 0.00999 0.04500 0.54000 0.58000 -17824.72077 0.53476 132 4 130 66 0.00159 \
0.00437 11.36141 2.12704831 -2617.99005 10.16739 10.16664 0.00407 0.00104 3417

Run TFA and BLS on the light curves listed in EXAMPLES/lc_list_tfa (the only light curve in that list is EXAMPLES/3.transit). For TFA we use the trendlist given in EXAMPLES/trendlist_tfa (the second and third columns are the x and y positions for each star), we have TFA only include trend stars that are more than 25 pixels from the EXAMPLES/3.transit (the coordinates for this star are given in EXAMPLES/lc_list_tfa). We pass the corrected light curve to the subsequent commands and do not output the TFA coefficients or the TFA model. We then run BLS on the light curve cleaned by TFA. We search over a fixed q range of 0.01 to 0.1 and a fixed period range of 0.5 to 5.0 days. We search 2000 frequency points with 200 bins per period. We add +7 to the UT time to get the local time (used to determine the fraction of delta-chi2 that comes from one night), and report the top two peaks in the BLS spectrum. We do not output the periodogram, but do output the bls model light curve in the third column of EXAMPLES/OUTDIR1/3.transit.bls.model. We subtract the best fit box-car transit before passing the light curve to the next command. We also compute the rms before running TFA and after running BLS. The best period is in the column BLS_Period_1_2 (the _2 just denotes that this is the third command) and is 2.12298216 days in this case (with a delta-chi2 of -25264.50212 given in the BLS_deltaChi2_1_2 column), the second best period is 1.35250338 days. The SN for the first peak (72.07052) is the signal-to-noise measured in the BLS spectrum, the more commonly used signal to pink-noise ratio for the transit is 14.99298. Note that there are 146 points in transit, 4 observed transits and the estimated red and white noises of the light curve (after subtracting the model transit) are 0.00140 and 0.00407 mag respectively. The depth of the first peak is 0.01162 mag and the q is 0.03. Note that TFA reduces the RMS from 0.00542 mag (column RMS_0) to 0.00471 mag (column TFA_RMS_1), and subtracting the transit reduces it further to 0.00407 mag (column RMS_3).


History:

July 2, 2008: Modified to allow input from stdin.

May 27, 2008: Modified the -BLS and -BLSFixPer commands to output the number of points in transit, the number of transits, the number of points before and after transit, the red and white noise of the signal after removing the bls model and the signal to pink noise ratio of the transit.

May 25, 2008: Fixed seg. fault in initialize_tfa routine. Added examples.

Mar 26, 2008: Fixed memory leak in -aov command.

Nov 26, 2007: Added -SYSREM and -binlc commands and the -jdtol option.

Nov 19, 2007: Added -TFA and -fluxtomag commands.

Nov 15, 2007: Added -MandelAgolTransit, -GetLSAmpThresh, and -aov_harm commands.