deflare - remove flares from light curves using lc_clean or lc_sigma_clip
deflare infile outfile method [nsigma] [plot] [save] [rateaxis] [pattern] [good_color] [exclude_color] [minlength] [mean] [stddev] [scale] [minfrac] [verbose]
The deflare routine is a Python script which allows the user to run either lc_clean() or lc_sigma_clip() on an input light curve file for filtering flares. These routines are called from the CIAO contributed script lightcurves.py (ahelp -l lightcurves), and can be run separately within a Python interpreter, such as chips, sherpa, or ipython. It is not necessary to load the deflare routine into a Python interpreter; the script performs this step automatically when it is run from the Unix command line, within the CIAO environment.
The lc_clean() routine removes flares (and periods of anomalously-low count rates) from light curves using the same technique as used on the ACIS backgroud files; it is based on the lc_clean program created by Maxim Markevitch. For more information, see the CIAO lc_clean ahelp file (ahelp -l lc_clean).
The lc_sigma_clip() routine - which used to be called analyze_ltcrv() prior to CIAO 4.1 - uses a simple iterative sigma-clipping routine to detect and remove flares from a light curve. For more information, see the CIAO lc_sigma_clip ahelp file (ahelp -l lc_sigma_clip).
The lc_clean algorithm was designed to work with the ACIS blank-sky datasets. The routine calculates - or is given - a mean rate from which it deduces a minimum and maximum valid count rate. Since the calculations assume Poisson statistics all the time bins are assumed to have the same width (i.e. exposure value). or further information, follow the documentation links in the CIAO thread "The ACIS "Blank-Sky" Background Files".
The 'mean', 'clip', 'sigma', 'scale', and 'minfrac' options control the behavior of the routine, and are represented in deflare by the parameters 'mean', 'nsigma', 'stddev', 'scale', and 'minfrac', respectively. If 'mean' is set then the value is taken as the initial guess, otherwise the level is calculated from the data, with outliers clipped using the 'nsigma' parameter (which defaults to 3.0). This is not an iterative clipping scheme - as used in lc_sigma_clip() - since the variance of the lightcurve is assumed to be
var = mean-rate / bin-width
(where sqrt(var) represents 'sigma' in lc_clean and 'stddev' in deflare). Once the mean level has been calculated, the limits are given by (if 'stddev' is changed from its default value of None).
mean - nsigma * stddev
mean + nsigma * stddev
otherwise they are
mean / scale
mean * scale
where scale defaults to 1.2. The minfrac option determines the minimum fraction of bins that must remain "good" for the routine to finish.
The lc_sigma_clip routine performs an iterative sigma-clipping algorithm, removing those points that fall outside the range
mean - nsigma
mean + nsigma
at each iteration until all data points are within this range. This algorithm is robust but not perfect; it can easily "overclean" a noisy lightcurve, or fail to reject any points, and so should not be used blindly.
The two options for lc_sigma_clip are 'sigma' and 'minlength', and are represented by deflare parameters 'nsigma' and 'minlength', respectively. The 'nsigma' argument determines the value used to clip the count rate data at each iteration, and defaults to 3.0. The 'minlength' option, which defaults to 3, is the minimum number of consecutive bins that must all lie within the final rate limits for that range to be considered good.
How do we use the output GTI file?
If the outfile argument is set, then a GTI file will be created listing those time periods which are considered to be "good" by the algorithm. This file can be used to filter an event file using the following syntax (we assume here that the event file to be filtered is evt2.fits and the output of the routine is called evt.gti):
unix% dmcopy "evt2.fits[@evt.gti]" evt2_filt.fits
A warning about low count-rate lightcurves
All the algorithms exclude bins which have a zero count rate, whether it is because the exposure time of the bin is zero or because no counts were detected in this bin. This is because the algorithms are not designed for use with low count-rate data. Therefore the routines should be used with extreme care when applied to data for which zero-rate bins are expected (for example a very-faint source).
A note on times for lc_sigma_clip
When run with verbose set to 1, durations may be displayed alongside time filters. These durations refer to the length of the interval and do not account for
- the dead-time correction of the observation (DTCOR header keyword);
- exposure variations within the bin due to existing time filters (i.e. the time filters already applied to the observation).
The sure way to find out what the final exposure time will be is to create a GTI file, use it to filter the event file, and then look at the EXPOSURE value of the filtered file.
unix% deflare lc.fits
The user is prompted to enter a flare-filtering method: 'clean' for lc_clean and 'sigma' for lc_sigma_clip. The default parameters for the chosen method are used to analyze the light curve in the file lc.fits.
unix% deflare lc.fits sigma outfile=clip.gti plot=yes save=lcsc_plot
Since the outfile parameter is set, the routine creates a GTI file (clip.gti) that represents the "good" times calculated by the lc_sigma_clip algorithm. After writing out the file, the filtered-out times will be displayed on the plot using red cross-hatched regions. We also save a PostScript plot of the light curve to the file lcsc_plot.ps.
The output GTI file can be used to filter an event file (say evt2.fits) like so:
unix% dmcopy "evt2.fits[@clip.gti]" evt2_filt.fits
unix% deflare lc.fits sigma outfile=out.gti nsigma=2.0 minlength=5 plot=yes
Create a GTI file, out.gti, using the lc_sigma_clip method, where the number of sigma used to determine the clipping of the count rate data per iteration is 2.0, and the minimum number of consecutive bins that must all lie within the final rate limits for that range to be considered good is 5. The resulting light curve and histogram plot showing the 'good' count rates and excluded time intervals are plotted to the screen.
unix% deflare lc.fits sigma outfile=out.gti plot=yes good_color=orange exclude_color=blue verbose=1
Create a GTI file using lc_clean with the default filtering parameters. Plot the light curve and histogram using orange for the 'good' data points and blue for the excluded times.
unix% deflare lc.fits clean outfile=out.gti pattern=none plot=yes
Create a GTI file using lc_clean and plot the light curve and histogram, but do not use regions to indicate the excluded times in the count-rate plot.
unix% deflare lc.fits clean outfile=out.gti nsigma=1.6 plot=no
Create a GTI file using the lc_clean method, where the number of sigma used to clip the data to calculate the mean for filtering is changed from the default 3.0 to 1.6.
unix% deflare lc.fits clean outfile= plot=yes save=lc_clean_plot
Save the plot of the light curve filtered with lc_clean to a PostScript file named "lc_clean_plot.ps", using default filtering parameters, without creating a GTI file. Note that since no GTI file is specified, the saved plot will NOT display the excluded regions, which are defined by deflare parameters 'good_color' and 'pattern'.
Detailed Parameter Descriptions
Parameter=infile (file required filetype=input)
The filename of the light curve to analyze.
The lc_clean and lc_sigma_clip routines are designed to work with light curves that were created using the CIAO dmextract tool, run with the opt parameter set to ltc1; see the "FORMAT OF LIGHT CURVES" section of this document for more details.
Parameter=outfile (file required filetype=output default=None)
The name of the GTI file to create.
If the outfile argument is set, then a GTI file will be created listing those time periods which are considered to be "good" by the algorithm. This file can be used to filter an event file.
Parameter=method (string required default=None)
The lightcurves routine to use for filtering the input light curve.
Either 'clean' for lc_clean or 'sigma' for lc_sigma_clip.
Parameter=nsigma (real not required default=3.0)
Number of sigma about the mean count rate to use to clip outliers.
This is the lc_clean 'clip' parameter and the lc_sigma_clip 'sigma' parameter. In lc_sigma_clip, 'sigma' determines which points are to be rejected at each iteration (reject if they are more than 'sigma' times the standard deviation away from the mean). In lc_clean, if the 'mean' parameter is set to the default 'None', then the mean count rate is calculted from the input data, using a single sigma clip of the data (so the 'clip' pararameter gives the number of standard deviations to retain).
Parameter=plot (boolean not required default=no)
Plot light curve and histograms of values?
If plot is 'yes' then plots of the light curve and histograms of the values will be displayed. If a GTI outfile has not been specified, then the plot will not display excluded time intervals (as they won't be calculated in this case). The plot may be saved to a PostScript file with the 'save' parameter, whether 'plot' is set to 'yes' or 'no'.
Parameter=save (file not required filetype=output default=None)
Name of PostScript file to which to save light curve and histogram plot.
The light curve and histogram plot produced by deflare (whether printed to the screen or not with the 'plot' parameter) will be saved to a PostScript file of the name specified in the 'save' parameter, unless set to the default "None".
Parameter=rateaxis (string not required default=y)
Should count-rate be on 'y' or 'x' axis of top plot?
The rateaxis option determines whether the light curve plot is drawn as count-rate on the y axis and time on the x axis (rateaxis='y') or the other way around (rateaxis='x').
Parameter=pattern (string not required default=solid)
Pattern to use in plot for filled regions representing excluded time intervals.
If a GTI outfile is specified and plot is 'True', then the time intervals excluded by the calculated GTI file are drawn on the plot as filled regions, using the pattern and color determined by the 'pattern' and 'pcol' options (setting pattern to 'none' stops this display).
The choices for this parameter are: nofill, solid, updiagonal, downdiagonal, horizontal, vertical, crisscross, brick, grid, hexagon, zigzag, wave.
Parameter=good_color (string not required default=green)
Color to use in plot to draw 'good' data points.
The 'good_color' parameter gives the color to use to draw points/histograms of the points that are considered good by the lc_clean and lc_sigma_clip algorithms (ignoring the minlength filter in lc_sigma_clip).
Parameter=exclude_color (string not required default=red)
Color to use in plot for excluded time intervals
If a GTI output file is specified and plot is 'True', then the time intervals excluded by the calculated GTI file are drawn on the plot as filled regions, using the 'pattern' and color determined by the 'pattern' and 'good_color' options.
Parameter=minlength (real not required default=3)
Minimum number of consecutive time bins which must pass rate filter in lc_sigma_clip routine.
The lc_sigma_clip minlength parameter determines the minimum number of consecutive time bins that must all pass the rate filter for the period to be considered 'good'. This parameter is not used within the sigma-clipping iteration. If minlength=1 then a filter on the count rate is determined instead.
Parameter=mean (real not required default=None)
Mean count rate (ct/s) of input light curve, as determined by lc_clean.
If 'mean' is 'None' (default), then the mean count rate is calculted from the input light curve, using a single sigma clip of the data (deflare 'nsigma' parameter, which is the same as lc_clean 'clip'). Otherwise, the mean count rate supplied should match the actual mean count rate of the input light curve. The standard deviation is calculated as sqrt(mean-rate/bin width). The mean value is used to clip the light curve.
Parameter=stddev (real not required default=None)
Standard deviation of signal, as determined by lc_clean.
The 'stddev' parameter is the lc_clean 'sigma' parameter, which is calculated from the 'mean' value as sqrt(mean-rate/bin width). When stddev is 'None', then the limits used for clipping the light curve are mean/scale to min*scale, otherwise they are mean - nsigma*stddev to mean + nsigma*stddev.
Parameter=scale (real not required default=1.2)
The scale factor about the mean count rate to use for light curve limits, as determined by lc_clean.
When the 'stddev' parameter is 'None', then the limits used by lc_clean for clipping the light curve are mean/scale to min*scale, otherwise they are mean - nsigma*stddev to mean + nsigma*stddev.
Parameter=minfrac (real not required default=0.1)
Minimum fraction of light curve bins that must pass the count-rate filter, as determined by lc_clean.
The 'minfrac' parameter gives the minimum fraction of the light curve bins that must pass the count-rate filter (after ignoring the zero count-rate bins).
Parameter=verbose (integer not required default=1 min=0 max=5)
Verbose can be from 0 to 5, generating different amounts of output.
Format of Lightcurves
The lc_clean routine is designed to work with light curves that were created using the CIAO dmextract tool, run with the opt parameter set to ltc1. They should however also work with any file which has the following columns:
- TIME_MIN and TIME_MAX (both columns are optional)
- COUNT_RATE or, if not present, RATE
If the file contains OBJECT and OBS_ID keywords then these will be used to label the plot.
Changes in the scripts 4.10.3 (October 2018) release
The ahelp file - the SYNTAX and PARAMETERS sections - has been updated to match the order of the parameters given in the deflare parameter file (e.g. 'plist deflare'). Please note that the exact set of parameters that are used depends on the method parameter.
Changes in the scripts 4.9.3 (May 2017) release
The script has been updated so that it will run under Python 3.5 (prior to this it could occasionally fail, depending on the options chosen and the input data).
Changes in the scripts 4.8.1 (December 2015) release
The code has been updated to avoid warning messages from NumPy version 1.9. There is no difference to how the script behaves.
About Contributed Software
This script is not an official part of the CIAO release but is made available as "contributed" software via the CIAO scripts page. Please see the installation instructions page for help on installing the package.
There are no known bugs for this tool.