diff --git a/docs/changelog.rst b/docs/changelog.rst new file mode 100644 index 0000000..ae43d27 --- /dev/null +++ b/docs/changelog.rst @@ -0,0 +1,18 @@ +=========== +Change Log +=========== + +v0.6 (2023-10-27) +--------------- + +* Now requires astropy regions >0.7 +* Lots of bug fixes +* Updated validate_det1_region for regions metadata change +* Added grouping keywords to spectral products +* Updates wrappers for stray light to add flexibility +* Added nustar_gen.diagnostic which has tools for determining if there is excess solar background and downloading GOES X-ray plots. +* Added additional example notebooks + +Contributors to this update: +- Brian Grefenstette +- Yash Bhargava diff --git a/docs/index.rst b/docs/index.rst index 52156c4..dc34b6d 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -5,6 +5,7 @@ nustar_gen documentation :maxdepth: 1 readme + changelog info diagnostic utils diff --git a/docs/io.rst b/docs/io.rst index 6c2c453..611771f 100644 --- a/docs/io.rst +++ b/docs/io.rst @@ -1,6 +1,6 @@ .. module:: nustar_gen -IO Functions +I/O Functions ================================== .. autofunction:: nustar_gen.io.find_be_arf diff --git a/docs/utils.rst b/docs/utils.rst index 60c2e33..fc4dbf6 100644 --- a/docs/utils.rst +++ b/docs/utils.rst @@ -7,4 +7,3 @@ Utility Scripts .. autofunction:: nustar_gen.utils.chan_to_energy .. autofunction:: nustar_gen.utils.make_usr_gti .. autofunction:: nustar_gen.utils.append_fits_entry - diff --git a/nustar_gen/diagnostic.py b/nustar_gen/diagnostic.py index 0923f1f..6a1c3dd 100644 --- a/nustar_gen/diagnostic.py +++ b/nustar_gen/diagnostic.py @@ -198,6 +198,7 @@ def goes_lightcurve(obs, show_sun=False, show_sky=False, show_impact=True): ax.add_patch(rect) plt.show() + return def compare_sun_spec(obs,mod='A',src_rad = 2*u.arcmin, en_range = (3, 20), en_bins = 34, ratio = False): @@ -221,7 +222,7 @@ def compare_sun_spec(obs,mod='A',src_rad = 2*u.arcmin, Energy range to use. Default is (3, 20). en_bins : int Number of bins to use for histogram (default is 34) - ratio : bool + ratio : bool Whether to plot the spectra or their ratio. Default is False. Returns @@ -313,7 +314,7 @@ def compare_sun_spec(obs,mod='A',src_rad = 2*u.arcmin, ax.axhline(1.0, color = 'green') plt.show() - + return \ No newline at end of file diff --git a/nustar_gen/io.py b/nustar_gen/io.py index 98c94b3..c12f709 100644 --- a/nustar_gen/io.py +++ b/nustar_gen/io.py @@ -12,7 +12,8 @@ def find_be_arf(): Returns ------- - Full path to the Be window file in the CALDB + be_arf_file : str + Full path to the Be window file in the CALDB Examples -------- @@ -48,8 +49,8 @@ def find_arf_template(): Returns ------- - Full path to the Be window file, which should be included when this repo is - checked out + arf_file : str + Full path to the Be window file, which should be included when this repo is checked out Examples -------- diff --git a/nustar_gen/utils.py b/nustar_gen/utils.py index de1f375..120daaa 100644 --- a/nustar_gen/utils.py +++ b/nustar_gen/utils.py @@ -243,7 +243,6 @@ def make_straylight_arf(det1im, regfile, filt_file, mod, obs): ''' Produces an ARF for a given observation. Currently uses counts-weighting to determine the contribution of the detabs parameters for each detector. - SHOULD PROBABLY BE MOVED TO WRAPPERS OR SOMETHING? Parameters ---------- @@ -345,6 +344,18 @@ def validate_det1_region(regfile): """ Code for making sure that region files that you're trying to use on the DET1 analysis code are in "image" coordinates + + Parameters + ---------- + regfile : str + Full path to a ds9 region file + + + + Returns + ------- + None + """ err=-1 import regions @@ -368,10 +379,11 @@ def validate_det1_region(regfile): # Check to make sure tha the first region in the file is an "include" region for ri in reg: - assert ri.meta['include'] is True, \ + assert ri.meta['include'] == 1, \ f'\n {regfile} has an exclusion region first! \n Put the source region first instead!' break + return def straylight_background(det1im_file='None', sky2det_file='None', reg_file = 'None', det1_expo_file = 'None', @@ -398,9 +410,6 @@ def straylight_background(det1im_file='None', sky2det_file='None', det1_expo_file : str Path to the DET1 effective area file. Should be the DET1 output of wrappers.make_exposure_map - - Optional Parameters - ------------------- diag : bool Show diagnostic information and make some plots showing what is happening. diff --git a/nustar_gen/version.py b/nustar_gen/version.py index 5103dc6..bd49665 100644 --- a/nustar_gen/version.py +++ b/nustar_gen/version.py @@ -5,4 +5,4 @@ from setuptools_scm import get_version version = get_version(root='..', relative_to=__file__) except Exception: - version = '0.6.dev9+g80e934a.d20220720' + version = '0.6.dev34+g013771b.d20231025' diff --git a/nustar_gen/wrappers.py b/nustar_gen/wrappers.py index ee07e9c..ef53165 100644 --- a/nustar_gen/wrappers.py +++ b/nustar_gen/wrappers.py @@ -4,15 +4,101 @@ from nustar_gen.utils import energy_to_chan, validate_det1_region from astropy import units as u +def make_image(infile, elow = 3, ehigh = 20, clobber=True, outpath=False, usrgti=False): + ''' + Spawn an xselect instance that produces the image in the energy range. + + Parameters + ---------- + infile: str + Full path to the file that you want to process + elow: float + Low-energy band for the image. Default is 3 keV. + ehigh: float + High-energy band for the image. Default is 20 keV. + + Other Parameters + ---------------- + + clobber: boolean, optional, default=True + Overwrite existing files? + + outpath: str, optional, default=os.path.dirname(infile) + Set the destination for output. Defaults to same location as infile. + + usrgti : str, optional, default = False + Use a GTI file to time-fitler the data (see nustar_gen.utils.make_usr_gti) + If False, do nothing. + + Return + ------- + outfile: str + The full path to the output image. + + ''' + # Make sure environment is set up properly + _check_environment() + + # Check if input file exists: + + # Account for the case where you only get the filename (i.e. if you're + # in the current working directory) + infile = os.path.abspath(infile) + + try: + with open(infile) as f: + pass + except IOError: + raise IOError("make_image: File does not exist %s" % (infile)) + + if not outpath: + outdir=os.path.dirname(infile) + else: + outdir=outpath + try: + os.makedirs(outdir) + except FileExistsError: + # directory already exists + pass + + # Trim the filename: + sname=os.path.basename(infile) + if sname.endswith('.gz'): + sname = os.path.splitext(sname)[0] + sname = os.path.splitext(sname)[0] + + + if usrgti is not False: + rshort = os.path.basename(usrgti) + rname = os.path.splitext(rshort)[0] + sname += f'_{rname}' + + + # Generate outfile name + outfile = os.path.join(outdir, sname+f'_{elow}to{ehigh}keV.fits') + logfile = os.path.join(outdir, sname+f'_{elow}to{ehigh}keV.log') + + if (os.path.exists(outfile)): + if clobber is True: + os.system("rm "+outfile) + else: + warnings.warn('make_image: %s exists, use clobber=True to regenerate' % (outfile)) + xsel_file = _make_xselect_commands(infile, outfile, elow, ehigh, usrgti=usrgti) + os.system(f"xselect @{xsel_file} > {logfile} 2>&1") +# os.system("rm -r -f "+xsel_file) + + return outfile + def make_spectra(infile, mod, src_reg, usrgti=False, mode='01', bgd_reg='None', outpath='None', runmkarf='yes', extended='no', grouping_flag=False, group_min_count = 30, group_pi_bad_high =1909, group_pi_bad_low =35, group_append='grp', oa_hist = False): ''' Generate a script to run nuproducts to extract a source (and optionally - a background) spectrum along with their response files. + a background) spectrum along with their response files. Does not automatically + run this script (must be run from the command line) separately. - Always runs numkrmf and numkarf for now. + Always runs numkrmf. Parameters ---------- @@ -30,21 +116,23 @@ def make_spectra(infile, mod, src_reg, usrgti=False, Other Parameters ------------------- - - bgd_reg: str + usrgti : str, optional, default=False + If not False must be the full path to the GTI that you want to apply + to the spectrum. + mode: str, optional, default = '01' + Optional. Used primarily if you're doing mode06 analysis and need to specify + output names that are more complicated. + bgd_reg: str, optional, default='None' If not 'None', then must be the full path to the background region file - - barycorr: bool - outpath: str - Optional. Default is to put the lightcurves in the same location as infile - - usrgti: str - Optional. Full path to the GTI that you want to apply to the spectrum. - - mode: str - Optional. Used primarily if you're doing mode06 analysis and need to specify - output names that are more complicated. + Output path for the resulting spectra. Default is to put the spectra in + the same location as infile. + runmkarf : str, optional, default='yes' + Flag for whether or not to generate an ARF. Default is 'yes'. + Alternative is 'no'. + extended : str, optional, default='no' + Flag for whether or not to use the extended ARF. Default is 'no'. + Alternative is 'yes'. grouping_flag: bool Optional. Used to indicate if the grouping of the spectra should be done or not group_min_count: int default value: 30 @@ -56,6 +144,10 @@ def make_spectra(infile, mod, src_reg, usrgti=False, group_append: str default value:'grp' A subtring which will be appended to the grouped spectra + Returns + -------- + scr : str + Path to the resulting shell script ''' @@ -145,7 +237,8 @@ def make_lightcurve(infile, mod, src_reg, barycorr=False, time_bin=100*u.s, mode='01',usrgti=False, bgd_reg='None', outpath='None', elow=3, ehigh=20): ''' - Generate a script to run nuproducts + Generate a script to run nuproducts to make a lightcurve. Script must be run + from the command line. Parameters ---------- @@ -163,28 +256,29 @@ def make_lightcurve(infile, mod, src_reg, Other Parameters ------------------- - bgd_reg: str - If not 'None', then must be the full path to the background region file - - barycorr: bool - Default is 'False'. If 'True', then queries the infile for the OBJ J2000 + barycorr: bool, optional, default=False + If True, then queries the infile for the OBJ J2000 coordinates and uses these for the barycenter correction. - - elow: float - Low-eneryg bound. Default is 3 keV. - - ehigh: float - High-energy bound. Default is 20 keV. - - time_bin: astropy unit - Length of the time bin. Default is 100*u.s. Min is 1*u.s - - outpath: str - Optional. Default is to put the lightcurves in the same location as infile - - mode: str - Optional. Used primarily if you're doing mode06 analysis and need to specify + time_bin: astropy unit, options, default=100 + Length of the time bin. Min is 1*u.s + mode: str, optional, default = '01' + Used primarily if you're doing mode06 analysis and need to specify output names that are more complicated. + usrgti : str, optional, default=False + If not False, must be the full path to the GTI that you want to apply to the spectrum. + bgd_reg: str, optional, Default='None' + If not 'None', then must be the full path to the background region file + outpath: str, optional + Optional. Default is to put the lightcurves in the same location as infile + elow: float, optional, default = 3 + Low-energy bound in keV + ehigh: float, optional, default = 20 + High-energy bound in keV + + Returns + -------- + scr : str + Path to the resulting shell script ''' @@ -297,6 +391,11 @@ def make_exposure_map(obs, mod, vign_energy = False, det_expo : boolean, optional, default=False Whether or not to retain the DET1 exposure map file + Returns + -------- + scr : str + Path to the script + ''' import glob @@ -355,96 +454,12 @@ def make_exposure_map(obs, mod, vign_energy = False, return expo_script -def make_image(infile, elow = 3, ehigh = 20, clobber=True, outpath=False, usrgti=False): - ''' - Spawn an xselect instance that produces the image in the energy range. - - Parameters - ---------- - infile: str - Full path to the file that you want to process - elow: float - Low-energy band for the image. Default is 3 keV. - ehigh: float - High-energy band for the image. Default is 20 keV. - - Other Parameters - ---------------- - - clobber: boolean, optional, default=True - Overwrite existing files? - - outpath: str, optional, default=os.path.dirname(infile) - Set the destination for output. Defaults to same location as infile. - - usrgti : str, optional, default = False - Use a GTI file to time-fitler the data (see nustar_gen.utils.make_usr_gti) - If False, do nothing. - - Return - ------- - outfile: str - The full path to the output image. - - ''' - # Make sure environment is set up properly - _check_environment() - - # Check if input file exists: - - # Account for the case where you only get the filename (i.e. if you're - # in the current working directory) - infile = os.path.abspath(infile) - - try: - with open(infile) as f: - pass - except IOError: - raise IOError("make_image: File does not exist %s" % (infile)) - - if not outpath: - outdir=os.path.dirname(infile) - else: - outdir=outpath - try: - os.makedirs(outdir) - except FileExistsError: - # directory already exists - pass - - # Trim the filename: - sname=os.path.basename(infile) - if sname.endswith('.gz'): - sname = os.path.splitext(sname)[0] - sname = os.path.splitext(sname)[0] - - - if usrgti is not False: - rshort = os.path.basename(usrgti) - rname = os.path.splitext(rshort)[0] - sname += f'_{rname}' - - - # Generate outfile name - outfile = os.path.join(outdir, sname+f'_{elow}to{ehigh}keV.fits') - logfile = os.path.join(outdir, sname+f'_{elow}to{ehigh}keV.log') - - if (os.path.exists(outfile)): - if clobber is True: - os.system("rm "+outfile) - else: - warnings.warn('make_image: %s exists, use clobber=True to regenerate' % (outfile)) - xsel_file = _make_xselect_commands(infile, outfile, elow, ehigh, usrgti=usrgti) - os.system(f"xselect @{xsel_file} > {logfile} 2>&1") -# os.system("rm -r -f "+xsel_file) - - return outfile def extract_sky_events(infile, regfile, elow=1.6, ehigh=165., clobber=True, outpath=False): ''' - Spawn an xselect instance that produces a new event file screened using a sky ds9 - region file. + Spawn an xselect instance that produces a new event file screened + using a sky ds9 region file. Parameters ---------- @@ -457,14 +472,12 @@ def extract_sky_events(infile, regfile, elow=1.6, ehigh=165., clobber=True, outp Other Parameters ---------------- - elow : float - Low energy cut to use. Default is 1.6 keV (PI == 0) - ehigh : float - High energy cut to use. Default is 165 keV. - + elow : float, optional, default=3 + Low energy cut to use. + ehigh : float, optional, default=20 + High energy cut to use. clobber: boolean, optional, default=True Overwrite existing files? - outpath: str, optional, default=os.path.dirname(infile) Set the destination for output. Defaults to same location as infile. @@ -546,9 +559,9 @@ def barycenter_events(obs, infile, mod='A', barycorr_pars=None): Other Parameters ------------------- - barycorr_pars: dict - additional parameters to be passed to barycorr, - in the format {'refframe': 'icrs', 'ra': 123.456} + barycorr_pars: dict, optional, default=None + If not None, must be a dict giving additional parameters to be + passed to barycorr, in the format {'refframe': 'icrs', 'ra': 123.456} ''' diff --git a/requirements.txt b/requirements.txt index b7b2f66..743978d 100644 --- a/requirements.txt +++ b/requirements.txt @@ -5,4 +5,4 @@ h5py scikit-image sphinx sphinx_rtd_theme -regions==0.5 +regions>=0.7