Reference for Mike functions
add_plot_unit
Synopsis
Create a new X-unit,Y-unit pair for isis_fancy_plots package
Usage
add_plot_unit(String_Type, String_Type [; ]);
Qualifiers
- SEE BELOW
Description
add_plot_unit(xunit, yunit; xscale=val, yscale=val, is_energy=val, xlabel=str, ylabel=str, pgfont=str);
Create a new X-unit and Y-unit pair for plot_counts, plot_data, plot_unfold, and plot_fit_model. (yunit and yscale only affect plot_unfold; plot_counts and plot_data are only sensitive to the choice of xunit.) Set is_energy=1 if the xunit scales as keV (otherwise it defaults to scaling as Angstroms). xscale is the scaling of keV/xunit or A/xunit. yscale scales *all* of the y-axes in the plot_unfold plots, so it must be set appropriately to achieve a desired effect. xunit and yunit will be treated as lower case regardless of input. Note that existing X-units will not be overwritten; their already defined values will be used regardless of input. A set of default axis labels will be produced, with "ylabel" substituting for:
"Photons cm^{-2} s^{-1}" (the default)
However, these labels can be rewritten using the new_plot_labels command.
To get the scalings correct, remember that plot_unfold(...;power=2) is always photons/s/cm^2 in the absence of scalings, power=0 is proportional to flux for wavelength based X-units, while power=3 is proportional to flux for energy based units. For these latter two cases, the yscale needs to account for divisions/multiplications of the scaled X-unit.
Examples:
To produce mJy vs. THz plots, use:
add_plot_unit("thz", "mjy"; xscale=Const_eV/Const_h*1.e12, yscale=Const_h*1.e26, xlabel="THz", ylabel="mJy/THz");
Const_eV/Const_h*1.e12 ~ 4.14e-3 is keV per Terra-Hz plot_unfold(...;power=2) produces photons/s/cm^2 => (photons*xunit)/s/cm^2/xunit (xunit scaling cancels out)
We wish to convert this to mJy:
(photons*xunit/s/cm^2/xunit) = [photons*(1.e-26 ergs)*Hz/(1.e-26 ergs)]/s/cm^2/Hz = Hz/(1.e-26 ergs)*mJy
hence multiplying by yscale=(1.e-26 ergs)/Hz = 1.e26*Const_h produces the desired scaling for plot_unfold.
To produce BTU/Acre/s vs. keV plots, use:
plot_unfold(...;power=3) for X-unit=keV produces keV/s/cm^2, which we wish to convert to BTU/Acre/s, hence xscale=1 and:
erg_p_btu = 1.055056e10 % ISO ergs/BTU cmsq_p_a = 4.0468564224e7 % cm^2 per international acre yscale = cmsq_p_a*Const_eV*1.e3/erg_p_btu % cm^2/Acre*BTU/keV
To achieve this same Y-unit vs. other energy based X-units, multiply yscale by keV/xunit, i.e., the value of xscale.
See also: fancy_plot_unit, set_plot_labels, new_plot_labels
apj_size
Synopsis
Set the pgplot output size to something suitable for ApJ single column (isis_fancy_plots package)
Usage
apj_size;
Description
Use as: isis> id = open_print("fig1.ps/vcps"); apj_size; nice_width; isis> plot(x,y); isis> close_print(id,"gv");
See also: sov, keynote_size, nice_width, open_print, close_print, pg_color, pg_info
append_chain
Synopsis
Append a chain fits file to another chain fits file
Usage
append_chain(String_Type chainfile1, String_Type chainfile2, String_Type chainoutfile)
Qualifiers
- verbose: show progress
Description
This function takes two chains stored in chainfile1 and chainfile2 written with the write_chain function (or by emcee itself) and appends the second one to the first one, again writing a fits file. This function does not make use of read_chain or write_chain, it's consisting only of fits routines. The function itself does perform some sanity and safety checks (same data, same fit_fun) but you have to make sure that you use the correct chains in the correct order.
See also: combine_chain, read_chain, write_chain, emcee
clear_all
Synopsis
Deletes data & corresponding ARFs and RMFs
Usage
clear_all();
Qualifiers
- noprompt avoid the prompt in scripts
Description
Deletes all data, the corresponding ARFs and RMFS, as well as the intrinsic correction factors used for Fermi spectra. This function removes the indicated spectra from the internal list; it does not affect the disk files containing the spectra.
See also: delete_data
close_print
Synopsis
Wrapper around close_plot to allow system function to call the output file (isis_fancy_plots package)
Usage
close_print(window_id, String_Type);
Description
Use as: isis> id = open_print("fig1.ps/vcps"); keynote_size; nice_width; isis> plot(x,y); isis> close_print(id,"gv");
See also: open_print, sov, open_plot, close_plot, apj_size, keynote_size, nice_width, pg_color, pg_info
combine_chain
Synopsis
Combine several mcmc chains generated with the same configuration
Usage
combine_chain(Array_Type chaininfiles, String_Type chainoutfile)
Qualifiers
- verbose: show progress
Description
This function takes mcmc chains stored in files passed as an array of strings in chaininfiles, written with the write_chain function (or by emcee itself) and combines them to one chain file, again writing a fits file. This function does not make use of read_chain or write_chain, it's consisting only of fits routines. The function itself does perform some sanity and safety checks (same data, same fit_fun) but you have to make sure that you use the correct chains in the correct order.
See also: read_chain, write_chain, emcee, append_chain
emcee_chain_hist_collect
Synopsis
Collect chain histograms from multiple emcee files
Usage
Struct_Type[] chist = emcee_chain_hist_collect( String_Type[] InFiles, Interger_Type[] PID )
or
```c
Struct_Type[] chist = emcee_chain_hist_collect( String_Type[] InFiles )
##### Qualifiers
* chatty: If given, enable informative output.
* ncut[=0]: #sim. steps to be cut from the start of the chains.
* autocut: If given, 'ncut' is automatically set, s.t. the faulty sim. steps
containing ZERO entries are deleted. These entries corresponds
to the initialisation of the walkers.
* id[=0]: Index of the InFiles to check consistency to.
* nbin[=100]: Number of histogram bins.
* pmin: Double_Type[length(PID)]: Containing parmater lower limits.
* pmax: Double_Type[length(PID)]: Containing parmater upper limits.
##### Description
Generate parameter histograms by collecting chains from several FITS Infiles
created with emcee. It is required that all InFiles are based on the same model
and data. Only InFiles are having the same 'model' key and same 'free_par'
indices are taken into account.
IMPORTANT is that a common random number server was used for all InFiles to
ensure there statistical independency! This is not checked automatically!
Chains from InFiles (passed the check) for each parameter in 'pid' are collected,
that is added to a common histogram with 'nbin' bins. With 'pmin' and 'pmax'
the limits of the histograms can be set manually, where 'pmin'/'pmax' must
be an array of the same length as 'pid'.
'id' sets the index in InFiles, which is used to perform the consistency check.
With 'ncut' the number of simulation iterations at the beginning of the
the chains which should be cut away are set.
'autocut' can be used to set 'ncut' automatically, s.t. all iteration
steps related to the walker initialization are cut away.
__See also__: emcee
#### emcee_merge
##### Synopsis
Merge emcee fits files
##### Usage
```c
emcee_merge( String_Type[] InFiles, String_Type OutFile)
or
```c
emcee_merge( String_Type[] InFiles )
##### Qualifiers
* chatty: If given, enable informative output.
* force: If given, overwrite existing OutFile.
* adjuststeps: If given, do not disregard InFiles with different #simsteps.
Instead cut all InFiles to min(#simsteps).
* ncut[=0]: #sim. steps to be cut from the start of the chains.
* autocut: If given, 'ncut' is automatically set, s.t. the faulty sim. steps
containing ZERO entries are deleted. These entries corresponds
to the initialisation of the walkers.
* id[=0]: Index of the InFiles to check consistency to.
##### Description
Merge several FITS Infiles created with emcee into one OutFile. It is required that
all InFiles are based on the same model and data. Only InFiles are merged that
have the same 'model' key and same 'free_par' indices compared to the first file
in InFiles.
IMPORTANT is that a common random number server was used for all InFiles to
ensure their statistical independency! This is not checked automatically!
InFiles (passed the check) are merged in a sense that the different walkers
from all InFiles are incorporated into a common chain, e.g., if InFiles are
two files with 5 free Parameters, 11 walkers and 1000 nsteps, respectively,
then the OutFile has 5 free Parameters, 22 walkers and 1000 nsteps.
To be able to merge this way NSTEPS in all InFiles has to be the same. InFiles
with different NSTEPS will be disregarded. With the 'adjuststeps' qualifier
it is possible to cut all InFiles to the mininmal occuring NSTEPS. NOTE that
this means to throw away simulation iterations in all InFiles.
With 'ncut' the number of simulation iterations at the beginning of the
the chains which should be cut away are set.
An existing OutFile can be overwritten with the 'force' qualifier.
If no OutFile is given, 'string_intersection' is used to determine the
OutFile path/name!
__See also__: emcee
#### fancy_plot_unit
##### Synopsis
Change the x-axis, and possibly the y-axis units, in the isis_fancy_plots package.
##### Usage
```c
fancy_plot_unit( String_Type [, String_Type]);
Description
fancy_plot_unit(xunit [,yunit]);
Change the X-axis plot units to "xunit" (as for the ISIS command plot_unit; and change the Y-axis unit to "yunit" (default yunit= "photons"). These will be used for the functions: plot_counts, plot_data, plot_unfold, plot_fit_model. Unit names are case insensitive.
Available X-units:
eV, keV, MeV, GeV, TeV, Angstrom, A, nm, um, mm, cm, m, Hz, kHz, MHz, GHz, psd (used for plotting power spectra from SITAR)
Available Y-units:
photons (default), mJy, ergs, watts, psd_leahy, psd_rms
Units added via add_plot_unit are also supported.
**NOTE**: Y-units photons/mJy/ergs/watts affect only plot_unfold, while psd_* are for plot_counts, but will also affect plot_data/plot_unfold.
Fundamentally, power=1 is proportional to photons/cm^2/s/xunit (plot_unfold) or Counts/bin (plot_counts), with higher (lower) powers multiplying (dividing) by xunit. plot_data is always Counts/sec/xunit ("power" has no effect).
mJy: y-unit for plot_unfold/power=2 is mJy.
ergs: y-unit for plot_unfold/power=3 (xunit=keV, etc.) or power=1 (xunit=A, etc.) is ergs/cm^2/sec.
watts: Similar behavior to the ergs unit, but yielding Watts/cm^2.
psd_leahy/psd_rms are for use with SITAR timing routines, and plot Power Spectra in Leahy or (RMS/Hz)^2 units vs. Hz, using plot_counts.
See also: plot_unit, add_plot_unit, set_plot_labels, new_plot_labels
keynote_size
Synopsis
Set the pgplot output size to something suitable for Keynoe presentations (isis_fancy_plots package)
Usage
apj_size;
Description
Use as: isis> id = open_print("fig1.ps/vcps"); keynote_size; nice_width; isis> plot(x,y); isis> close_print(id,"gv");
See also: sov, apj_size, nice_width, open_print, close_print, pg_color, pg_info
load_scaled_data
Synopsis
Loads a dataset with a non-unity AREASCAL keyword or column
Usage
id = load_scaled_data( String_Type file);
Qualifiers
- rmf [="rmf.fits", response file if not in data header]
- arf [="arf.fits", effective area file if not in data header]
- bkg [="bkg.fits", background file if not in data header]
Description
Load a data file with a non-unity AREASCAL, rewrite the data/background to undo the ISIS default application of this value, create a post_model_hook to properly apply it, and rewrite the data and background backscales to properly apply it. The data itself must either reference a set of responses/backgrounds, or these must be entered via qualifiers upon initial loading of the data.
See also: load_data, _define_back, set_data_backscale
new_plot_labels
Synopsis
Create new plot labels for routines in isis_fancy_plots package.
Usage
new_plot_labels(String_Type [, String_Type]);
Qualifiers
- xlabel : New X-axis label
- ylabel : String array with new Y-axis labels for data
- rlabel : String array with new Y-axis labels for residuals
- clabel : String array with new Y-axis residual labels for Cash statistics
- mllabel : String array with new Y-axis residual labels for Maximum Likelihood statistics
- vlabel : String with new Doppler velocity X-axis label
- zlabel : String with new redshift X-axis label
- pg_font : PGPLOT font type (default is \fr)
Description
new_plot_labels(x_unit [,y_unit];xlabel=string,ylabel=[string],..., pg_font=string);
Changes the default labels for plot_counts, plot_data, plot_unfold for the different units set by fancy_plot_unit. pg_font string will be *prepended* to all inputs. Use any valid pair of units to set residual, Cash statistic, or velocity/redshift axis labels.
Use:
set_plot_labels(;pg_font="\fr") -or- set_plot_labels(;pg_font=\fr
)
to restore defaults.
Inputs:
x_unit : angstrom,a,nm,um,micron,cm,mm,m, ev,kev,mev,gev,tev, hz,khz,mhz,ghz, psd y_unit : photons (=default), ergs, watts, mjy, psd_leahy, or psd_rms xlabel : String with new X-axis label ylabel : String *array* (up to 9 elements) with new Y-axis labels. Order of the array *must* be: [plot_data, plot_unfold(power=0->3),plot_counts(power=0->3)] rlabel : String *array* (up to 3 elements) with new Y-axis labels for residuals. Order of the the array *must* be: [chi, chi^2, ratio] clabel : String Array with new Cash statistic labels [res=1 or 4, 2 or 5, which yield +/-sqrt(|Delta C|), Delta C, respectively] mllabel: String Array with new Maximum Likelihood statistic labels [res=1 or 4, 2 or 5, which yield +/-sqrt(|Delta ML|), Delta ML, respectively] vlabel : String with new Doppler velocity X-axis label zlabel : String with new Redshift X-axis label pg_font: "\fn", "\fr", "\fi", "\fs" = normal, roman, italic, or script fonts will be used on the labels (Default is \fr.)
See also: set_plot_labels, fancy_plot_unit, add_plot_unit
nice_width
Synopsis
Sets reasonable defaults for plot line widths.
Usage
nice_width;
Description
Equivalent to: isis> set_plot_widths(;d_width=2, de_width=2, r_width=2, re_width=2, m_width=2);
See also: set_plot_widths, apj_size, keynote_size, open_print, close_print
open_print
Synopsis
Wrapper around open_plot (and pg_color) to allow a system function to call the output file (isis_fancy_plots package)
Usage
id = open_print(String_Type);
Description
Use as: isis> id = open_print("fig1.ps/vcps"); keynote_size; nice_width; isis> plot(x,y); isis> close_print(id,"gv");
See also: close_print, sov, open_plot, close_plot, apj_size, keynote_size, nice_width, pg_color, pg_info
pg_color
Synopsis
Set pgplot colors 17, 18, 19, 20 to a green, brown, pink, dark yellow (isis_fancy_plots package)
Usage
pg_color; -or- pgcolor;
Description
See also: pg_info, pginfo
pg_info
Synopsis
Print a core dump of some useful pgplot and isis_fancy_plots information.
Usage
pg_info; -or- pginfo;
Qualifiers
- NONE
Description
See also: Nearly all isis_fancy_plot functions return a use message if invoked without arguments
plot_comps
Synopsis
Create a data plot with model components explicitly shown (isis_fancy_plots package)
Description
plot_comps({data},pstruct,&plot_func); % pstruct=struct{dcol, mcol, ...} plot_comps({data},&plot_func;dcol={val},mcol={val},ccol={val},cstyle=val,...);
Use a fancy plotting function, e.g., plot_counts or plot_data or plot unfold, passed as a reference, and cycle through all the components with a norm parameter. Plot each of these as a separate model component. The plot functions now take two additional optional qualifiers (which alternatively can be passed via the pstruct structure variable): ccol and cstyle. The ccol parameter gives the color of the model components for each dataset, which can be different from the color of the model for the complete model. The cstyle allows a global change of the line_style for *all* of the model components. (I.e., only one alternate line_style can be chosen.) data is the usual combination of integers (=individual data sets), arrays (=data sets to be combined), and lists (=id of combined datasets). Examples:
plot_comps({1,[2,3]},popt,&plot_counts;xrange={1,10}); plot_comps({5,8},&plot_unfold); plot_comps{{{1}},popt,&plot_data);
See also: plot_counts, plot_data, plot_unfold, plot_residuals, plot_fit_model, plotxy, plot_comps, plot_double
plot_counts
Synopsis
Plot counts per bin (isis_fancy_plots package)
Description
plot_counts({indx,[arry],{cid}},pstruct); % pstruct = struct{ dcol, mcol, rcol, ...} plot_counts({indx,[arry],{cid}};dcol={val},mcol={val},rcol={val,[arry],val},...);
Plot background subtracted data, model, and residuals as counts/bin Residuals are units of chi, chi2, or ratio, and will be based upon whether one chooses sigma=model, data, or gehrels in set_fit_statistic(); (data error bars are only affected by the latter two). set_fit_statistic("cash"); will alter the residuals to the Cash statistic. set_fit_statistic("ml"); will alter the residuals to the Maximum Likelihood statistic.
Options below refer to structure variables/qualifiers
indx = list of data set indices to be plotted. Any indices grouped in an array within that list will be *combined* in the data plot. Single number in list is combo id, {#} = [combination_members(#)]. dcol = (pgplot) color value for data (or list of color values) decol = (pgplot) color value for data error bars (or list of color values) mcol = (pgplot) color value for model (or list of color values) 0 => No model plotted rcol = (pgplot) color value for residuals (or list of color values; arrays within the list allow for individual color values if portions of the data are combined, but their associated residuals are not)
recol = color for residual error bars (or list of color values; arrays within the list act as for residual color inputs) dsym = (pgplot) symbol value for data (or list of symbol values) 0 => histogram plot rsym = (pgplot) symbol value for residuals (or list of symbol values; arrays within the list act as for residual color inputs) 0 => histogram plot xrange = List of X-limits for the data & model & residuals Note: Any X- or Y-range set to NULL is autoscaled yrange = List of Y-limits for the data & model and (optionally) residuals oplt = 0 (default) for new plot, !=0 for overplotting no_reset= 0 (default)- plots *will* be reset, i.e., next plot moves to new pane (multiplot), next plot redraws window (single plot). no_reset=1 necessary for overplotting multiplots (oplt=1 sufficient for single plots). res = 0 (default), no residuals; 1, 2, or 3 = chi, chi2, or ratio residuals 4, 5, or 6 = chi, chi2, or ratio, but combine residuals for combined data set_fit_method("cash"); or set_fit_method("ml") will cause res=(2 or 5) or res=(1 or 4) to plot the residual for the Cash or Maximum Likelihood statistic or its square root, respectively power = 0, 1 (default), 2, or 3 for Counts/bin X (1/Unit, 1, Unit, Unit^2), respectively bkg = List of 0's (subtract background-default), 1's (include backgrounds), or -1's (plot *only* the background [no model plotted in this case]). Ratio residuals will include background in data/model, other residuals are unaffected. Indices within a combination are treated the same.) xlabel = String that will overwrite default X-axis label (default=NULL) ylabel = String or string array that will overwrite the default Y-axis labels (second element of array applies to residuals; default=NULL) zshift = List of redshifts to be applied to the data (default zshift={0,0,...}) vzero = If set, the reference X-unit value to be defined as zero velocity. The X-axis then becomes a velocity axis (km/s) referenced to this point (default vzero=NULL; setting vzero/zaxis supersedes zshift) zaxis = If not 0, use a redshift axis instead of a velocity axis *if* vzero is defined (default zaxis=0) scale = Multiplicatively scale the Y-axis by the values in a list. Any arrays in the list should hold the individual scalings for data set arrays in the input index list. **Note:** these values only scale the plots, not the fits. (Default values are 1.) gap = 1 (default), models are histograms with gaps where data has gaps, 0 , model are bin-centered lines, without gaps.
See also: plot_data, plot_unfold, plot_residuals, plot_fit_model, plotxy, plot_comps, plot_double
plot_data
Synopsis
Plot counts per unit per second (isis_fancy_plots package)
Description
plot_data({indx,[arry],{cid}},pstruct); % pstruct = struct{ dcol, mcol, rcol, ...} plot_data({indx,[arry],{cid}};dcol={val},mcol={val},rcol={val,[arry],val},...);
Plot background subtracted data, model, and residuals as counts/xunit/sec. Residuals are units of chi, chi2, or ratio, and will be based upon whether one chooses sigma=model, data, or gehrels in set_fit_statistic(); (data error bars are only affected by the latter two). set_fit_statistic("cash"); will alter the residuals to the Cash statistic. set_fit_statistic("ml"); will alter the residuals to the Maximum Likelihood statistic.
Options below refer to structure variables/qualifiers
indx = list of data set indices to be plotted. Any indices grouped in an array within that list will be *combined* in the data plot. Single number in list is combo id, {#} = [combination_members(#)]. dcol = (pgplot) color value for data (or list of color values) decol = (pgplot) color value for data error bars (or list of color values) mcol = (pgplot) color value for model (or list of color values) 0 => No model plotted rcol = (pgplot) color value for residuals (or list of color values; arrays within the list allow for individual color values if portions of the data are combined, but their associated residuals are not) recol = color for residual error bars (or list of color values; arrays within the list act as for residual color inputs) dsym = (pgplot) symbol value for data (or list of symbol values) 0 => histogram plot rsym = (pgplot) symbol value for residuals (or list of symbol values; arrays within the list act as for residual color inputs) 0 => histogram plot xrange = List of X-limits for the data & model & residuals Note: Any X- or Y-range set to NULL is autoscaled yrange = List of Y-limits for the data & model and (optionally) residuals oplt = 0 (default) for new plot, !=0 for overplotting no_reset= 0 (default)- plots *will* be reset, i.e., next plot moves to new pane (multiplot), next plot redraws window (single plot). no_reset=1 necessary for overplotting multiplots (oplt=1 sufficient for single plots). res = 0 (default), no residuals; 1, 2, or 3 = chi, chi2, or ratio residuals 4, 5, or 6 = chi, chi2, or ratio, but combine residuals for combined data set_fit_method("cash"); or set_fit_method("ml") will cause res=(2 or 5) or res=(1 or 4) to plot the residual for the Cash or Maximum Likelihood statistic or its square root, respectively bkg = List of 0's (subtract background-default), 1's (include backgrounds), or -1's (plot *only* the background [no model plotted in this case]). Ratio residuals will include background in data/model, other residuals are unaffected. Indices within a combination are treated the same.) xlabel = String that will overwrite default X-axis label (default=NULL) ylabel = String or string array that will overwrite the default Y-axis labels (second element of array applies to residuals; default=NULL) zshift = List of redshifts to be applied to the data (default zshift={0,0,...}) vzero = If set, the reference X-unit value to be defined as zero velocity. The X-axis then becomes a velocity axis (km/s) referenced to this point (default vzero=NULL; setting vzero/zaxis supersedes zshift) zaxis = If not 0, use a redshift axis instead of a velocity axis *if* vzero is defined (default zaxis=0) scale = Multiplicatively scale the Y-axis by the values in a list. Any arrays in the list should hold the individual scalings for data set arrays in the input index list. **Note:** these values only scale the plots, not the fits. (Default values are 1.) sum_exp = If==1, then when combining data sets, sum the exposure times (as opposed to using the mean exposure time; default sum_exp=1). gap = 1 (default), models are histograms with gaps where data has gaps, 0 , model are bin-centered lines, without gaps.
See also: plot_data, plot_unfold, plot_residuals, plot_fit_model, plotxy, plot_comps, plot_double
plot_double
Synopsis
Use two different plot functions in the same figure (isis_fancy_plots package)
Description
plot_double({data},pstruct,&plot_funcI,&plot_funcII); % pstruct=struct{dcol, mcol, ...} plot_double({data},&plot_funcI,&plot_funcII;dcol={val},mcol={val},ccol={val},...);
Using fancy plotting functions, e.g., plot_counts or plot_data or plot unfold, passed as references, apply plot_funcI in the upper panel and plot_funcII in the lower panel. (If residuals are chosen to be displayed, include a third panel for them beneath the first two plots.) In each plot, cycle through all the components with a norm parameter. Plot each of these as a separate model component. The plot functions now take two additional optional qualifiers (which alternatively can be passed via the pstruct structure variable): ccol and cstyle. The ccol parameter gives the color of the model components for each dataset, which can be different from the color of the model for the complete model. The cstyle allows a global change of the line_style for *all* of the model components. (I.e., only one alternate line_style can be chosen.) data is the usual combination of integers (=individual data sets), arrays (=data sets to be combined), and lists (=id of combined datasets). Plot parameters retain their usual meaning, with the exception of yrange and power. yrange is now a list of up to 6 elements, with the first two applying to plot_funcI, the next two applying to plot_funcII, and the final two applying to the residuals. If power is a list of two elements, then the first element applies to plot_funcI and the second element applies to plot_funcII. Examples:
plot_double({1,[2,3]},popt,&plot_unfold,&plot_counts;xrange={1,10}); plot_double({5,8},&plot_unfold,&plot_data); plot_double{{{1}},popt,&plot_data,&plot_counts);
See also: plot_counts, plot_data, plot_unfold, plot_residuals, plot_fit_model, plotxy, plot_comps, plot_double
plot_fit_model
Synopsis
Plot background subtracted model as counts/xunit/sec (isis_fancy_plots package).
Description
plot_fit_model({indx,[arry],{cid}},pstruct); % pstruct = struct{ mcol, ...} plot_fit_model({indx,[arry],{cid}};mcol={val},...);
Plot background subtracted model as counts/xunit/sec.
Options below refer to structure variables/qualifiers
indx = list of data set indices to be plotted. Any indices grouped in an array within that list will be *combined* in the data plot. Single number in list is combo id, {#} = [combination_members(#)]. mcol = (pgplot) color value for model (or list of color values) 0 => No model plotted xrange = List of X-limits for the data & model & residuals Note: Any X- or Y-range set to NULL is autoscaled yrange = List of Y-limits for the data & model and (optionally) residuals oplt = 0 (default) for new plot, !=0 for overplotting no_reset= 0 (default)- plots *will* be reset, i.e., next plot moves to new pane (multiplot), next plot redraws window (single plot). no_reset=1 necessary for overplotting multiplots (oplt=1 sufficient for single plots). bkg = List of 0's (subtract background-default), 1's (include backgrounds), or -1's (plot *only* the background [no model plotted in this case]). Ratio residuals will include background in data/model, other residuals are unaffected. Indices within a combination are treated the same.) xlabel = String that will overwrite default X-axis label (default=NULL) ylabel = String or string array that will overwrite the default Y-axis labels (second element of array applies to residuals; default=NULL) zshift = List of redshifts to be applied to the data (default zshift={0,0,...}) vzero = If set, the reference X-unit value to be defined as zero velocity. The X-axis then becomes a velocity axis (km/s) referenced to this point (default vzero=NULL; setting vzero/zaxis supersedes zshift) zaxis = If not 0, use a redshift axis instead of a velocity axis *if* vzero is defined (default zaxis=0) scale = Multiplicatively scale the Y-axis by the values in a list. Any arrays in the list should hold the individual scalings for data set arrays in the input index list. **Note:** these values only scale the plots, not the fits. (Default values are 1.) sum_exp = If==1, then when combining data sets, sum the exposure times (as opposed to using the mean exposure time; default sum_exp=1). gap = 1 (default), models are histograms with gaps where data has gaps, 0 , model are bin-centered lines, without gaps.
See also: plot_counts; plot_data, plot_unfold, plot_residuals, plotxy, plot_comps, plot_double
plot_residuals
Synopsis
Plot the data residuals (isis_fancy_plots package)
Description
plot_residuals({indx,[arry],{cid}},pstruct); % pstruct = struct{ rcol, rsym, ...} plot_residuals({indx,[arry],{cid}};dcol={val},rcol={val,[arry],val},rsym=...);
Plot data residuals, without plotting the data. Residuals related to data indices in [arry] will appear in a single ascii file if write_plot(); is used, and will be combined in the plot if res>3 is chosen.
Options below refer to structure variables/qualifiers
indx = list of data set indices to be plotted. Any indices grouped in an array within that list will be *combined* in the data plot. Single number in list is combo id, {#} = [combination_members(#)]. rcol = (pgplot) color value for residuals (or list of color values; arrays within the list allow for individual color values if portions of the data are combined, but their associated residuals are not) recol = color for residual error bars (or list of color values; arrays within the list act as for residual color inputs) rsym = (pgplot) symbol value for residuals (or list of symbol values; arrays within the list act as for residual color inputs) 0 => histogram plot xrange = List of X-limits for the data & model & residuals Note: Any X- or Y-range set to NULL is autoscaled yrange = List of Y-limits for the data & model and (optionally) residuals oplt = 0 (default) for new plot, !=0 for overplotting no_reset= 0 (default)- plots *will* be reset, i.e., next plot moves to new pane (multiplot), next plot redraws window (single plot). no_reset=1 necessary for overplotting multiplots (oplt=1 sufficient for single plots). res = 0 (default), no residuals; 1, 2, or 3 = chi, chi2, or ratio residuals 4, 5, or 6 = chi, chi2, or ratio, but combine residuals for combined data set_fit_method("cash"); or set_fit_method("ml") will cause res=(2 or 5) or res=(1 or 4) to plot the residual for the Cash or Maximum Likelihood statistic or its square root, respectively xlabel = String that will overwrite default X-axis label (default=NULL) ylabel = String or string array that will overwrite the default Y-axis labels (second element of array applies to residuals; default=NULL) zshift = List of redshifts to be applied to the data (default zshift={0,0,...}) vzero = If set, the reference X-unit value to be defined as zero velocity. The X-axis then becomes a velocity axis (km/s) referenced to this point (default vzero=NULL; setting vzero/zaxis supersedes zshift) zaxis = If not 0, use a redshift axis instead of a velocity axis *if* vzero is defined (default zaxis=0)
See also: plot_counts, plot_data, plot_unfold, plot_fit_model, plotxy, plot_comps, plot_double
plot_unfold
Synopsis
Plot flux-corrected spectra (isis_fancy_plots package).
Description
plot_unfold({indx,[arry],{cid}},pstruct); % pstruct = struct{ dcol, mcol, rcol, ...} plot_unfold({indx,[arry],{cid}};dcol={val},mcol={val},rcol={val,[arry],val},...);
Plot background subtracted unfolded data, model, and residuals using a variety of X- and Y-units set by fancy_plot_unit(xunit, [yunit]); Residuals are units of chi, chi2, or ratio, and will be based upon whether one chooses sigma=model, data, or gehrels in set_fit_statistic(); (data error bars are only affected by the latter two). set_fit_statistic("cash"); will alter the residuals to the Cash statistic. set_fit_statistic("ml"); will alter the residuals to the Maximum Likelihood statistic.
Options below refer to structure variables/qualifiers
indx = list of data set indices to be plotted. Any indices grouped in an array within that list will be *combined* in the data plot. Single number in list is combo id, {#} = [combination_members(#)]. dcol = (pgplot) color value for data (or list of color values) decol = (pgplot) color value for data error bars (or list of color values) mcol = (pgplot) color value for model (or list of color values) 0 => No model plotted rcol = (pgplot) color value for residuals (or list of color values; arrays within the list allow for individual color values if portions of the data are combined, but their associated residuals are not) recol = color for residual error bars (or list of color values; arrays within the list act as for residual color inputs) dsym = (pgplot) symbol value for data (or list of symbol values) 0 => histogram plot rsym = (pgplot) symbol value for residuals (or list of symbol values; arrays within the list act as for residual color inputs) 0 => histogram plot xrange = List of X-limits for the data & model & residuals Note: Any X- or Y-range set to NULL is autoscaled yrange = List of Y-limits for the data & model and (optionally) residuals oplt = 0 (default) for new plot, !=0 for overplotting no_reset= 0 (default)- plots *will* be reset, i.e., next plot moves to new pane (multiplot), next plot redraws window (single plot). no_reset=1 necessary for overplotting multiplots (oplt=1 sufficient for single plots). res = 0 (default), no residuals; 1, 2, or 3 = chi, chi2, or ratio residuals 4, 5, or 6 = chi, chi2, or ratio, but combine residuals for combined data set_fit_method("cash"); or set_fit_method("ml") will cause res=(2 or 5) or res=(1 or 4) to plot the residual for the Cash or Maximum Likelihood statistic or its square root, respectively power = 0, 1 (usual default), 2 (default for mJy) or 3 (default for ergs/Watts vs. energy units)=> photons/cm^2/s/xunit *(1/xunit,1,xunit,xunit^2) bkg = List of 0's (subtract background-default), 1's (include backgrounds), or -1's (plot *only* the background [no model plotted in this case]). Ratio residuals will include background in data/model, other residuals are unaffected. Indices within a combination are treated the same.) xlabel = String that will overwrite default X-axis label (default=NULL) ylabel = String or string array that will overwrite the default Y-axis labels (second element of array applies to residuals; default=NULL) zshift = List of redshifts to be applied to the data (default zshift={0,0,...}) vzero = If set, the reference X-unit value to be defined as zero velocity. The X-axis then becomes a velocity axis (km/s) referenced to this point (default vzero=NULL; setting vzero/zaxis supersedes zshift) zaxis = If not 0, use a redshift axis instead of a velocity axis *if* vzero is defined (default zaxis=0) scale = Multiplicatively scale the Y-axis by the values in a list. Any arrays in the list should hold the individual scalings for data set arrays in the input index list. **Note:** these values only scale the plots, not the fits. (Default values are 1.) gap = 1 (default), models are histograms with gaps where data has gaps, 0 , model are bin-centered lines, without gaps. con_mod= 1 (default), the smear the model by the detector response, otherwise plot the unsmeared model at the internal resolution of the arf
Note: Model flux is: ( \int dE S(E) )/dE, while data is (C(h) - B(h))/(\int R(h,E) A(E) dE)/dh/t, where A(E) is effective area, R(h,E) is RMF, C(h)/B(h) are total/background counts. Thus, the data and model will match best only in the limit of a delta function RMF, and in fact might look different than the residuals (which is the only proper comparison between data and model, anyhow).
See also: set_power_scale, plot_counts, plot_data, plot_residuals, plot_fit_model, plotxy, plot_comps, plot_double
plotxy
Synopsis
Generate a simple x/y plot with error bars (isis_fancy_plots package)
Description
plotxy(x,dxm,dxp,y,dym,dyp, pstruct); % pstruct = struct{dcol, decol, xrng, ...} plotxy(x,dxm,dxp,y,dym,dyp; dcol=#, decol=#, ...);
Also accepts:
plotxy(x,y [,pstruct;qualifiers]); plotxy(x,dxm,dxp,y [,pstruct;qualifiers]); plotxy(x,,,y,dym,dyp [,pstruct;qualifiers]);
Plot simple x,y plots with error bars: x-dxm, x+dxp, etc.
Options below refer to structure variables/associative keys/qualifiers:
dcol = (pgplot) color value for data decol = (pgplot) color value for data dsym = (pgplot) symbol value for data Note: dsym=0 will *not* produce a histogram plot xrange= List of X-limits for the data. Ranges previously input via xrange(); will be respected if this option is not set. yrange= List of Y-limits for the data. Ranges previously input via yrange(); will be respected if this is not set. xlabel= String for the X-axis label ylabel= String for the Y-axis label Note: xlabel(); ylabel(); commands will also work. oplt = 0 (default) for new plot, !=0 for overplotting
Further note: plotxy(...); will apply the choices from connect_points(#);
See also: plot_counts, plot_data, plot_unfold, plot_residuals, plot_fit_model, plotxy, plot_comps, plot_double
reset_plot_defaults
Synopsis
Changes some of the defaults on the isis_fancy_plots package
Usage
reset_plot_defaults(;dcol=Integer_Type, ...);
Qualifiers
- dcol : Data color value
- dsym : Data symbol value
- mcol : Model color value
- sum_exp : !=0, Sum the exposures when combining data
- use_con_flux : !=0, the unfolded model includes response smearing
- gap : ==0, plot models across data gaps.
Description
reset_plot_defaults(; dcol=#, dsym=#, mcol=#, sum_exp=#, use_con_flux=#, gap=#);
Resets some of the plot defaults to a user's specifications. (See help messages for individual plotting functions, and pg_info, to understand these settings.) Use with no arguments to see current values.
See also: pg_info
set_plot_labels
Synopsis
Restore default plot labels of isis_fancy_plots package
Usage
set_plot_labels(;pg_font="\\fr") -or- set_plot_labels(;pg_font=``\fr``)
Qualifiers
- pg_font : ="\fn", "\fr", "\fs", or "\fi"
Description
See also: new_plot_labels, fancy_plot_unit, add_plot_unit
set_plot_widths
Synopsis
Sets plot widths for isis_fancy_plots package
Usage
set_plot_widths([;qualifiers]);
Qualifiers
- m_width : Model line width
- d_width : Data line width
- de_width : Data error bar line width
- r_width : Residual line width
- re_width : Residual error bar line width
- ebar_x : X error bar term cap length
- ebar_y : Y error bar term cap length
- data_err : !=0 X error bars plotted
Description
set_plot_widths(; d_width=#, de_width=#, r_width=#, re_width=#, m_width=#, ebar_x=#, ebar_y=#, data_err=#);
Sets line widths on data, residuals, error bars, and models, sets the length of x/y error bar term caps (ebar_x, ebar_y), and toggles the X error bar plotting. Values are retained until explicitly overwritten.
See also: nice_width, pg_info
set_power_scale
Synopsis
Determine the y-axis is scaled by energy when using the isis_fancy_plot package.
Usage
set_power_scale(Integer_Type);
Description
set_power_scale(a);
Determine how the y-axis is scaled by energy or wavelength when using the fancy_plot routines. a=1, E|lambda is set to the midpoint of the bin. a=2, E|lambda is set to the geometric midpoint (i.e., sqrt(Elo*Ehi)). a=3, E|lambda is set to the average value assuming an x^-2 powerlaw. Any other values, and this message is printed.
See also: fancy_plot_unit, add_plot_unit, plot_unfold, plot_counts
sitar_avg_cpd
Synopsis
Create an averaged cross power spectral density.
Usage
(f,psd_a,psd_b,cpd,navg,avg_a,avg_b,) = sitar_avg_cpd(cnts_a,cnts_b,l);
Qualifiers
- dt [Length of evenly spaced bins. Default == 1.]
- times [Times of measurements (otherwise presumed to have no gaps).]
- norm [Determine PSD normalization convention. =1, 'Leahy normalization'. Poisson noise PSD == 2, in absence of deadtime effects. !=1, rms or Belloni-Hasinger or Miyamoto normalization, i.e., PSD == (rms)^2/Hz & Poisson noise== 2/Rate.]
Description
Take an evenly spaced lightcurve (presumed counts vs. time), and calculate the PSD and CPD in segments of length l, averaged over the whole lightcurve. Segments with data gaps are skipped. Deadtime corrections to Poisson noise are *not* made. Poisson noise is *not* subtracted from average PSDs.
Inputs: cnts_a/b: Arrays of total counts in each time bin l : Length of individual PSD segments (use a power of 2!!!)
Outputs: f : PSD frequencies ( == 1/Input Time Unit) psd_a/b : Average PSDs cpd : Normalized Cross Power Spectral Density (complex array) navg : Number of data segments going into the averages avg_a/b : Average number of counts per segment of length l
See also: sitar_avg_psd, sitar_lbin_cpd, sitar_lbin_psd, sitar_define_psd, sitar_lags
sitar_avg_psd
Synopsis
Create an averaged power spectral density.
Usage
(f,psd,navg,avg_cnts [,psd_err]) = sitar_avg_psd(cnts,l);
Qualifiers
- dt [Length of evenly spaced bins. Default == 1.]
- times [Times of measurements (otherwise presumed to have no gaps).]
- norm [Determine PSD normalization convention. =1, Leahy normalization. Poisson noise PSD == 2, in absence of deadtime effects. !=1, rms or Belloni-Hasinger or Miyamoto normalization, i.e., PSD == (rms)^2/Hz & Poisson noise== 2/Rate.]
- err [If it exists, also return the PSD error calculated directly by assuming each sample has 100% error.]
Description
Take an evenly spaced lightcurve (presumed counts vs. time), and calculate the PSD in segments of length l, averaged over the whole lightcurve. Segments with data gaps are skipped. Deadtime corrections to Poisson noise are *not* made. Poisson noise is *not* subtracted from average PSD.
Inputs: cnts : Array of total counts in each time bin l : Length of individual PSD segments (use a power of 2!!!)
Outputs: f : PSD frequencies ( == 1/Input Time Unit) psd : Average PSD. navg : Number of data segments going into the average avg_cnts: Average number of counts per segment of length l
Optional Output: psd_err : Estimated error on the PSD, calculated from lightcurve sample
See also: sitar_avg_cpd, sitar_lbin_cpd, sitar_lbin_psd, sitar_define_psd, sitar_lags
sitar_bin_events
Synopsis
Bin an event-based lightcurve.
Usage
lc = sitar_bin_events(t,dt,gti_lo,gti_hi);
Qualifiers
- tstart [Beginning of first output time bin. Default: min(t).]
- tstop [End of last output time bin. Default: max(t).]
- obin [For a tstop, last bin width might be < dt. Include this overflow bin if width is > obin*dt. Default: obin=0.1.]
- user_bin.lo [Lower time bounds for user defined grid.]
- user_bin.hi [Upper time bounds for a user defined grid. Overrides dt, tstart, and tstop inputs.]
- minexp [Do not include bins with exposure < minexp. Default: 0.]
- delgap [If set, delete 0 exposure bins from the lightcurve, otherwise set their rate & error to 0.]
Description
Inputs: t : Discrete times at which rates are measured dt : Width of evenly spaced time bins in rebinned lightcurve (Ignored if a user_bin is input.)
Outputs: lc.rate : Mean rate of resulting binning lc.err : Error on bin rate (minimum = -log(0.32)/(bin exposure)) lc.cts : Counts in a bin lc.bin_lo : Lower time bounds of binned lightcurve lc.bin_hi : Upper time bounds of binned lightcurve lc.expos : Exposure of a time bin (i.e., intersection with good time intervals)
See also: sitar_rebin_rate
sitar_define_psd
Synopsis
Store a power spectral density as a fittable dataset.
Usage
id = sitar_define_psd(flo,fhi,psd,err [,noise]);
Qualifiers
- noise [Same as optional input, noise.]
Description
Inputs: f_lo : Low value of PSD frequency bin (Hz -> keV) f_hi : High value of PSD frequency bin (Hz -> keV) psd : Array of PSD values (Power -> Counts/bin) err : Array of PSD Errors
Optional Input: noise : Array (*or single value*) of Poisson noise levels. If undefined, background is undefined.
Outputs: id : Dataset Index
See also: sitar_avg_cpd, sitar_lbin_cpd, sitar_avg_psd, sitar_define_psd, sitar_lags
sitar_global_optimum
Synopsis
Loads a dataset with a non-unity AREASCAL keyword or column
Usage
ans = sitar_global_optimum( cell, ncp_prior, type [; first, mean_prior, rate_prior, alpha=Double_Type, beta=Double_Type] );
Description
Inputs: cell : A structure, with cell.pops, cell.size, cell.lo_t cell.hi_t,cell.dtcor (see sitar_make_data_cells) ncp_prior : Parameter for prior on number of 'blocks' type : Identical to type from sitar_make_data_cells
Qualifiers: first : If set, the code runs in 'trigger' mode, and returns upon the first sign of a change, so long as it is greater than 'first' cells from the beginning of the lightcurve
EVENT MODE PRIORS (type=1 or 2)-
Default prior is that p_1, the probability of one or more photons in a frame, is uniformly dsitributed from 0->1.
alpha : >=0 implies p_1 prior is (1+alpha) (1-p_1)^alpha
mean_prior: if set, prior on *rate* is exp(-Lambda/Lambda_0), where Lambda_0 is the mean rate from the entire observation. Supercedes any choice on alpha.
rate_prior: if set, prior is chosen to be uniform between rate_low, rates ranging from rate_low -- rate_high, with rate_high defaults of rate_low = 1/3 mean rate and rate_high = 3 times mean rate. Supercedes any choice of alpha or mean_prior.
BINNED MODE PRIORS AND POSTERIORS (type=3)-
alpha, : The prior on the bin rate (sort of) goes as beta (Lambda)^(alpha-1) Exp(-beta*Lambda), where Lambda=rate, and the default is alpha=beta=1
max_like : Use maximum likelihood, log(Pmax) = N log(N/T) -N, as the posterior test function. Overrides any choice of alpha & beta for binned mode.
Outputs: results.cpt : Array of change point locations for the maximum likelihood solution (indices are specific to the cell input); results.last,: (Diagnostic purposes only) Arrays of the .best location of the last change point and the associated maximum log probability results.cts : Counts in each block results.rate : Rate in each block results.err : Poisson error for the block rate results.lo_t,: Times of lower (>=) and upper (<) block .hi_t boundaries.
See also: sitar_make_data_cells
sitar_glvary
Synopsis
Create an optimal lightcurve using the Gregory-Loredo algorithm.
Usage
gl = sitar_glvary(t);
Qualifiers
- tmin [The minimum time to consider. Default=min(t)-1]
- tmax [The maximum time to consider. Default=max(t)+1]
- mmax [Consider lightcurve divisions from 2 to mmax evenly spaced bins. Default=300.]
- thresh [Truncate the number of partitions by ignoring those for which \sum_m(odds ratio)<max(\sum_n(odds ratio))/exp(thresh), where n = 2 -> mmax. thresh sets a minimum probability: (1-p_min) ~ exp(thresh)*(1-p_peak). Default=2.]
- nbins [Create an output lightcurve with nbins. Default=mmax.]
- texp, frac_exp [Array pair that give for fractional exposure as a function of time, which will be interpolated and used to correct lightcurve rates. Arrays must have a minimum of five entries. Default is for no correction.]
Description
Use the Gregory-Loredo algorithm to find the odds ratios that even divisions of a lightcurve are better descriptions of the data than a constant lightcurve. A 'best estimate' lightcurve can also be output for a lightcurve with nbins.
Inputs: t : Array of event times
Outputs: gl.p : Total probability that some evenly partitioned lightcurve, with up to gl.mmax bins, is a better description than a constant lightcurve gl.ppart : The probability for an individual evenly partitioned lightcurve that it is a better description than a constant lightcurve gl.lodds_sum : The natural logarithm of the sum of the odds ratios comparing lightcurves with two or more partitions to a constant lightcurve. gl.p == exp(gl.lodds_sum)/[1+exp(gl.lodds_sum)] gl.mpeak : The number of partitions for the evenly partitioned lightcurve with the maximum probability. gl.mmax : The maximum number of partitions actually used (influenced by the setting of the thresh parameter) gl.m : The number of partitions corresponding to each evenly partitioned lightcurve considered (=[2:mmax]) gl.pm : Total probability that some evenly partitioned lightcurve is a better description than a constant lightcurve for each maximum number of partitions considered ([2:mmax]) gl.nj : The counts histogram corresponding to each partitioning above gl.aj : The integrated fractional exposure for each partitioning above gl.a_avg : The averaged fractional exposure. gl.tmin : The value of tmin actually used (maximum of [tmin,min(texp)]); gl.tmax : The value of tmax actually used (minimum of [tmax,max(texp)]); gl.tlc : The output lightcurve times (an array with input nbins bins) gl.rate : Best estimate of the lightcuve rates at the above times gl.erate : Best estimate of the lightcuve rate errors at the above times
See also: sitar_global_optimum
sitar_lags
Synopsis
Create time lags and coherence function from input power spectra.
Usage
(lag,dlag,g,dg) = sitar_lags(freq,psda,psdb,cpd,noisea,noiseb [,navg]);
Description
Given two input PSD (without noise subtracted), their CPD, and their associated noise levels, all as functions of Fourier frequency, calculate the time lag and coherence function (and associated errors) vs. f
Inputs: freq : Fourier frequency array psda/b : Power spectra array associated with two lightcurves (Noise subtraction *not* applied!) cpd : Cross power spectra array associated with two lightcurves noisea/b: Power spectra noise level. Array *or* a constant.
Optional Unput: navg : Number of averages (over independent data segments and frequency bins) associated with each input frequency bin. Can be a constant vs. f (defaulted to 1).
Outputs: lag : Time lag vs. frequency. Negative values mean psdb lags psda dlag : Associated error on the lag. g : Coherence function dg : Error on the coherence function
See also: sitar_avg_cpd, sitar_lbin_cpd, sitar_avg_psd, sitar_lbin_psd, sitar_define_psd
sitar_lbin_cpd
Synopsis
Logarithmically rebin a cross power spectral density.
Usage
(aflo,afhi,apsda,apsdb,acpd,nf,[anoisea,anoiseb]) = sitar_lbin_cpd(f,psda,psdb,cpd,dff,[noisea,noiseb,&rev]);
Description
Logarithmically rebin two PSDs (and, optionally, their associated Poisson noise levels) and the associated Cross Power Spectral Density
Inputs: f : Array of Fourier frequencies psda/b : Arrays of PSD values dff : Delta f/f value for logarithmic bin spacing
Optional Inputs: noisea/b : Arrays (*or single values*) of Poisson noise levels rev : If declared (i.e., 'isis> variable rev;') and input, rev returns the reverse indices for the binning (i.e., (apsda[i] = mean( psda[ rev[i] ] ), etc.)
Outputs: aflo : Lower boundary of Fourier frequency bin afhi : Upper boundary of Fourier frequency bin apsda/b : Rebinned Power Spectral Density values acpd : Rebinned Cross Power Spectral Density values nf : Number of frequencies going into the bin anoisea/b: Rebinned noise level (array, even for single value input)
See also: sitar_avg_cpd, sitar_avg_psd, sitar_lbin_psd, sitar_define_psd, sitar_lags
sitar_lbin_psd
Synopsis
Logarithmically rebin a power spectral density.
Usage
(aflo,afhi,apsd,nf [,anoise,apsd_err]) = sitar_lbin_psd(f,psd,dff);
Qualifiers
- noise [Array (*or single value*) of Poisson noise levels.]
- rev [If declared (i.e., 'isis> variable rev;') and input, rev returns the reverse indices for the binning (i.e., af[i] = mean( freq[ rev[i] ] ), etc.).]
- psd_err [Array of errors for the PSD, in which case a PSD error array will be returned.]
Description
Logarithmically rebin a PSD (and, optionally, its associated Poisson noise and error values)
Inputs: f : Array of Fourier frequencies psd : Array of PSD values dff : Delta f/f value for logarithmic bin spacing
Outputs: aflo : Lower bounds of Fourier frequency bins afhi : Upper bounds of Fourier frequency bins apsd : Rebinned PSD values nf : Number of frequencies going into the bin
Optional Outputs: anoise : Rebinned noise level (array, even for single value input) apsd_err: Rebinned PSD error array
See also: sitar_avg_cpd, sitar_lbin_cpd, sitar_avg_psd, sitar_define_psd, sitar_lags
sitar_make_data_cells
Synopsis
Loads a dataset with a non-unity AREASCAL keyword or column
Usage
cell = sitar_make_data_cells(tt,type,max_delt,frame,tstart,tstop);
Description
Create data cells from event data to be used with the Bayesian blocks algorithm. tt are the times of the events; type=1,2, or 3 determines cell breaks (midway between events; right before subsequent event, or binned); events closer than max_delt are grouped together in a single cell; output cell sizes are in units of frame (or = bin size for type=3), tstart/tstop are the stop/start times for the output cells.
Output is a structure with cell.pops (number of events per cell), cell.size (duration of cell in units for frame), cell.lo/hi_t (start/stop times of cell), cell.dtcor (array, currently set to unity, for storing cell deadtime corrections).
See also: sitar_global_optimum
sitar_pfold_event
Synopsis
Fold an event-based lightcurve on a given period (with derivatives).
Usage
profile = sitar_pfold_event(t, p, gti_lo, gti_hi);
Qualifiers
- pdot [Period derivative.]
- pddot [Period second derivative.]
- nphs [Number of phase bins in output folded lightcurve.]
- tzero [Time of zero phase. Default = t[0].]
- phase_lo [Arrays for phase bin boundaries (e.g., for uneven phase bins).]
- phase_hi [Arrays must be same length with matched boundaries. Supercedes nphs.]
Description
Inputs: t : Times at which rates are measured rate : Lightcurve rate p : Period on which to fold the lightcurve gti_lo : Array of lower boundaries of good time intervals gti_hi : Array of upper boundaries of good time intervals
Outputs: profile.bin_lo: Start value of phase bin profile.bin_hi: Stop value of phase bin profile.cts : Counts in phase bin profile.var : Square root of counts in phase bin profile.expos : Integrated exposure in phase bin
See also: sitar_pfold_rate, sitar_epfold_rate
sitar_pfold_rate
Synopsis
Fold a rate-based lightcurve on a given period (with derivatives).
Usage
profile = sitar_pfold_rate(t,rate,p);
Qualifiers
- pdot [Period derivative.]
- pddot [Period second derivative.]
- nphs [Number of phase bins in output folded lightcurve.]
- tzero [Time of zero phase. Default = t[0].]
- phase_lo [Arrays for phase bin boundaries (e.g., for uneven phase bins).]
- phase_hi [Arrays must be same length with matched boundaries. Supercedes nphs.]
Description
Inputs: t : Times at which rates are measured rate : Lightcurve rate p : Period on which to fold the lightcurve
Outputs: profile.bin_lo: Start value of phase bin profile.bin_hi: Stop value of phase bin profile.mean : Mean rate in phase bin profile.var : Variance of rate in phase bin profile.sdm : Standard deviation of mean rate in a phase bin profile.num : Number of data points in phase bin
See also: sitar_pfold_event, sitar_epfold_rate
sitar_readasm
Synopsis
Read and RXTE ASM data file
Usage
event = sitar_readasm( file );
Qualifiers
- tstart [Start time to be read; units of MET or MJD]
- tstop [Stop time to be read; units of MET or MJD]
- maxchi2 [Maximum acceptable chi2 for ASM solution]
- chnl [0 for total band; !=0 for three ASM colors+total]
- mjd [Changes from default of Mission Elapsed Time (days) to MJD]
- jd [Changes from default of Mission Elapsed Time (days) to JD]
Description
Inputs: file : ASM filename, e.g. 'xa_cygx1_d1.lc' or 'xa_cygx1_d1.col'
Outputs: event.time: Time of each event event.rate: Total ASM count-rate event.err : Uncertainty in count rate event.chi2: Chi2 of ASM solution
Optional Outputs: event.[ch1,ch2,ch3]_rate: ASM count rate in each channel event.[ch1,ch2,ch3]_err : ASM count rate error in each channel event.[ch1,ch2,ch3]_chi2: Chi2 for ASM solution in each channel (Overrides event.chi2, which won't be output)
sitar_rebin_rate
Synopsis
Rebin a rate lightcurve
Usage
lc = sitar_rebin_rate(t,dt,rate [,err]);
Qualifiers
- tstart [Beginning of first output time bin]
- tstop [End of last output time bin]
- minbin [Minimum required events per bin, else the bin is set to 0, or ignored. Default=1.]
- user_bin.lo [Lower time bounds for user defined grid]
- user_bin.hi [Upper time bounds for a user defined grid. Overrides dt, tstart, and tstop inputs]
- delgap [If set, delete empty bins from the lightcurve, otherwise, set them to 0]
- weight [If set and err input, the mean is weighted by the error. I.e., mean = (\sum rate/err^2)/(\sum 1/err^2), and the output variance becomes (\sum (mean-rate/err^2)^2)/(\sum 1/err^2)^2]
Description
Variables in [] are optional, but are order specific unless a qualifier.
Inputs: t : Discrete times at which rates are measured dt : Width of evenly spaced time bins in rebinned lightcurve (Ignored if a user_bin is input.) rate : (Presumed GTI & Exposure Corrected) rates
Optional Input: err : Error on the rates
Outputs: lc.rate : Mean rate of resulting rebinning lc.bin_lo : Lower time bounds of rebinned lightcurve lc.bin_hi : Upper time bounds of rebinned lightcurve lc.num : Number of events going into a bin lc.var : Variance of the rate in a time bin. (Only calculated where lc.num >= 2, otherwise set to zero.)
Optional Outputs: lc.err : Rate errors combined in quadrature, i.e., it's sqrt{\sum{ err^2 }}/N, where N is the number of points in that particular bin (i.e., lc.num). Only computed if err is input [otherwise, just use sqrt(lc.var)].
See also: sitar_bin_events
sov
Synopsis
Set pgplot outer viewport (isis_fancy_plots package)
Usage
sov(Double_Type, Double_Type, Double_Type, Double_Type);
Description
sov(xmin,xmax,ymin,ymax);
Equivalent to: isis> v=struct{xmin,xmax,ymin,ymax}; isis> v.xmin=xmin; isis> v.xmax=xmax; isis> v.ymin=ymin; isis> v.ymax=ymax; isis> set_outer_viewport(v) and: isis> set_outer_viewport(struct{xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax});
See also: set_outer_viewport, apj_size, keynote_size, nice_width, open_print, close_print, pg_color, pg_info
write_plot
Synopsis
Write the data from the last plot made with the isis_fancy_plots package
Usage
[pd =] write_plot(root; head=0, data);
Qualifiers
- head : ==0 to suppress header in output file
- data : if set, return data as a structure instead
Description
[pd =] write_plot(root; head=0, data);
Creates ASCII files with the data from the last plot made (with any plot device) using plot_counts, plot_data, plot_residuals, or plot_unfold. Data will be stored in columns in the following order- bin_lo, bin_hi, data_value, data_error (single column), model, {mean residual, residual-1sigma, residual+1sigma}. There can be multiple sets of columns for the groups of residuals if the data are combined, but the residuals are not.
Important notes:
The X-axis and Y-axis units will reflect those of the plot. The binning and noticed data ranges will reflect data filters applied within ISIS, but not the (usually more limited) plot range filters.
Each uncombined data set, or group of combined data sets, will be be written to a separate ASCII file. A file for combined data sets can have multiple columns for residuals if they are left uncombined.
Combined data or residuals will only write out the chosen combinations, not the individual pieces used to create them.
Files created from plot_residuals will only output the bin_lo/hi grid and the residual/-/+ columns, not the data or model.
If the model or residuals were not plotted, they will not be written to the ASCII file.
Plots made from plot_unfold will create separate files for the data and the unfolded model, as they can be on different grids. The model file will begin with the appropriate bin_lo, bin_hi.
Inputs:
root- A string with the base file name. Outputs will be-- root_#.dat, root_#.res (for plot residuals), and root_#.mod (for plot_unfold, if the data and model are on separate grids). # will will cycle from 0 to the number of input data groups-1, assigning data sets in the order that they were input to the plot command. ***Old files will be overwritten if 'root' is not a unique name.***
head- Set as a qualifier only. If !=0 (default=1), the ASCII data files will have useful header information appended, otherwise, the files will just contain columns of data.
data- If the data qualifier exists, instead of writing an ASCII file, the write_plot function will return a structure with the data values from the last plot generated.
See also: pg_info