Merge remote-tracking branch 'origin/main' into ics_sample

.github/workflows/ci_tests.yml

@@ -19,8 +19,8 @@ jobs:
     strategy:
       matrix:
         os: [ubuntu-latest]
-        python: ['3.10', '3.11']
-        deps: [current, numpy124, astropydev, numpydev, astropydev-numpydev]
+        python: ['3.11', '3.12']
+        deps: [current, numpy211, astropydev, numpydev, astropydev-numpydev]
 
     steps:
     - name: Check out repository
@@ -32,10 +32,10 @@ jobs:
     - name: Install base dependencies
       run: |
         python -m pip install --upgrade pip
-    - name: Test with numpy = 1.24
-      if: "contains(matrix.deps, 'numpy124')"
+    - name: Test with numpy = 2.1.1
+      if: "contains(matrix.deps, 'numpy211')"
       run: |
-        python -m pip install numpy==1.24
+        python -m pip install numpy==2.1.1
     - name: Test with dev version of numpy
       if: "contains(matrix.deps, 'numpydev')"
       run: |

bin/frb_pzdm_mag

@@ -0,0 +1,16 @@
+#!/usr/bin/env python
+#
+# See top-level LICENSE file for Copyright information
+#
+# -*- coding: utf-8 -*-
+
+
+"""
+This script performs estimates P(z|DM) and limiting mag
+"""
+
+from frb.scripts import pzdm_mag
+
+if __name__ == '__main__':
+    args = pzdm_mag.parser()
+    pzdm_mag.main(args)

docs/installing.rst

@@ -26,13 +26,13 @@ We recommend that you use `Anaconda <https://www.continuum.io/downloads/>`_
 to install and/or update these packages.
 
 * `python <http://www.python.org/>`_ versions 3.8 or later
-* `numpy <http://www.numpy.org/>`_ version 1.22 or later
+* `numpy <http://www.numpy.org/>`_ version 2.1 or later
 * `astropy <http://www.astropy.org/>`_ version 5.1 or later
 * `scipy <http://www.scipy.org/>`_ version 1.11 or later
 * `healpy <https://healpy.readthedocs.io/en/latest/index.html>`_ version 1.15 or later
 * `pandas <https://pandas.pydata.org/>`_ version 1.5 or later
 * `requests <https://pillow.readthedocs.io/en/5.3.x/>`_  version 2.18 or later
-* `extinction <https://extinction.readthedocs.io/en/latest/>`_ version 0.4.2 or greater
+* `extinction <https://extinction.readthedocs.io/en/latest/>`_ version 0.4.6 or greater
 * `matplotlib <https://matplotlib.org/>`_ version 3.7 or greater
 * `linetools <https://github.com/linetools/linetools>`_  version 0.3 or later
 * `astropath <https://github.com/FRBs/astropath>`_  version 0.1 or later
@@ -58,6 +58,7 @@ The following package(s) is/are required to access FRB galaxy spectra:
 The following package is required to map a slit onto a finder chart (frb.figures.finder):
 
 * `photutils <https://photutils.readthedocs.io/en/stable/>`_  version 1.11.0 or later
+* `scikit-image <https://scikit-image.org/>`_  version 0.21.0 or later
 
 The following are required to use our KCWI datacube handling tools:
 

docs/scripts.rst

@@ -107,11 +107,18 @@ This script takes as input the FRB DM and its coordinates (approximate
 are fine) and then estimates the redshift range assuming 
 the Macquart relation (and simple host + MW contributions, optionally 
 input).  For an (optionally input; tuple) confidence interval, 
-it reports back the putative redshift range for the FRB.
+it reports back the putative redshift range for the FRB. It also 
+allows for plotting the host redshift range on the magnitude vs redshift 
+evolution and setting a title for the figure.  These calculations can be 
+done assuming a few different telescope models (CHIME, DSA, Parkes, FAST,
+CRAFT, CRAFT_ICS_892/1300/1632) or a perfect telescope model (default). 
+The telescope models are used to determine the DM-z grids that have been 
+computed with the zdm code/repository.
 
 Here is the usage::
 
-    usage: frb_pz_dm [-h] [--dm_hostmw DM_HOSTMW] [--cl CL] coord DM_FRB
+    usage: frb_pz_dm [-h] [--dm_hostmw DM_HOSTMW] [--cl CL] coord DM_FRB 
+                    [--magdm_plot] [--fig_title FIG_TITLE] [--telescope TELESCOPE]
 
     Script to print a summary of an FRB to the screen [v1.0]
 
@@ -129,6 +136,16 @@ Here is the usage::
     --cl CL               Confidence limits for the z estimate [default is a 95
                             percent c.l., (2.5,97.5)]
 
+    --magdm_plot          Plot the host redshift range given DM on the magnitude
+                            vs redshift evolution. Default=False.   
+
+    --fig_title FIG_TITLE  title for the figure; e.g., FRBXXXXX
+
+    --telescope TELESCOPE  telescope model for the DM-z grid: CHIME, DSA, Parkes,
+                            FAST, CRAFT, CRAFT_ICS_892/1300/1632, perfect. Default
+                            = perfect
+
+
 
 frb_sightline
 =============

frb/dm/prob_dmz.py

@@ -1,7 +1,9 @@
 """ Code for calculations of P(DM|z) and P(z|DM)"""
 import numpy as np
 import os
-from pkg_resources import resource_filename
+
+
+from importlib.resources import files
 
 from scipy.stats import norm, lognorm
 from scipy.interpolate import interp1d, interp2d
@@ -150,32 +152,35 @@ def build_grid_for_repo(outfile:str):
     print("This will be used going forth")
 
 
-def grab_repo_grid():
+
+def grab_repo_grid(grid_name):
     """
-    Grab the grid from the Repository
-    This may require the code to build it first!
+    Grab the grid from the Repository based on the given grid name
+
+    Args:
+        grid_name (str): Name of the grid to grab
 
     Returns:
-        dict: Numpy dict from the npz save file
+        dict: Numpy dict from the npz or npy save file
     """
 
     # File
-    PDM_z_grid_file = os.path.join(
-        resource_filename('frb', 'data'), 'DM',
-        'PDM_z.npz')
-
+    grid_file = files('frb.data.DM').joinpath( grid_name)
+
     # Build?
-    if not os.path.isfile(PDM_z_grid_file):
-        build_grid_for_repo(PDM_z_grid_file)
-
+    if grid_name == 'PDM_z.npz':
+        if not os.path.isfile(grid_file):
+            build_grid_for_repo(grid_file)
+
     # Load
-    print(f"Loading P(DM,z) grid from {PDM_z_grid_file}")
-    sdict = np.load(PDM_z_grid_file)
+    print(f"Loading P(DM,z) grid from {grid_file}")
+    sdict = np.load(grid_file)
 
     # Return
     return sdict
 
 
+
 def get_DMcosmic_from_z(zarray, perc=0.5, redo_pdmz_grid=False, DMevals=np.linspace(1.,2000.,1000), beta=3., F=0.31, cosmo=defs.frb_cosmo):
     """
     Gives DMcosmic values of zarray, considering the percentile.

frb/galaxies/eazy.py

@@ -3,8 +3,8 @@
 import os
 import warnings
 from pkg_resources import resource_filename
-from distutils import spawn
 import subprocess
+import shutil
 
 import numpy as np
 import pandas
@@ -373,7 +373,7 @@ def run_eazy(input_dir, name, logfile):
     _, param_file, translate_file = eazy_filenames(input_dir, name)
 
     # Find the eazy executable
-    path_to_eazy = spawn.find_executable('eazy')
+    path_to_eazy = shutil.which('eazy')
     if path_to_eazy is None:
         raise ValueError("You must have eazy in your Unix path..")
     # Run it!

frb/galaxies/mag_dm.py

@@ -20,7 +20,7 @@
 
 
 def r_vs_dm_figure(z_min, z_max, z, PzDM, outfile='fig_r_vs_z.png',
-               flipy=True, known_hosts = False):
+               flipy=True, known_hosts = False, title=None, logz_scale=False):
     """
     Plots the intersection of galaxy apparent magnitude evolution with redshift 
     and redshift distribution of the FRB and saves it. 
@@ -40,6 +40,10 @@ def r_vs_dm_figure(z_min, z_max, z, PzDM, outfile='fig_r_vs_z.png',
     f_mL = frb_gal_u.load_f_mL()
 
 
+
+    # Increase default font size
+    plt.rcParams.update({'font.size': 12})
+
     # set up the figure
     plt.figure(figsize=(6, 5))
     gs = gridspec.GridSpec(1,1)
@@ -75,7 +79,9 @@ def r_vs_dm_figure(z_min, z_max, z, PzDM, outfile='fig_r_vs_z.png',
 
     # plot
     # L curves
-    zvals = np.linspace(0.021, 4.0, 200)
+
+
+    zvals = np.linspace(0.021, z_max+0.5, 200)
     m_Lstar = f_mL(zvals)
     ax.plot(zvals, m_Lstar, '-r', label='L*')
 
@@ -87,13 +93,14 @@ def r_vs_dm_figure(z_min, z_max, z, PzDM, outfile='fig_r_vs_z.png',
 
 
     # Add P(z|DM)
-    xmnx = (0., 4.)
-    ymnx = (14, 33.)
+    xmnx = (0., z_max+0.5)
+    ymnx = (14, max(15, m_001Lstar[-1]+0.5))
     if flipy:
         ymnx = (ymnx[1], ymnx[0])
 
     y_range = np.linspace(ymnx[0], ymnx[1], 300) 
     x_range = np.linspace(z_min,z_max,300)
+
     pzd = np.interp(x_range,z,PzDM)
     X,Y = np.meshgrid(x_range,y_range)
     Z = X/X*pzd
@@ -102,20 +109,28 @@ def r_vs_dm_figure(z_min, z_max, z, PzDM, outfile='fig_r_vs_z.png',
     #zlbl = 'P(z|DM) [95% c.l.]'
     c=ax.pcolor(X, Y, Z*1000, cmap='Blues')
     zlbl = '95% c.l. FRB Redshift \n estimated from DM'
-    ax.text(np.mean([z_min,z_max]), 20.5, zlbl,
+    text_x = 0.3 * (xmnx[1] - xmnx[0]) + xmnx[0]
+    text_y = 0.5 * (ymnx[1] - ymnx[0]) + ymnx[0]
+    ax.text(text_x, text_y, zlbl,
             color='k', 
             size='large', ha='center')
 
-    ax.set_xlabel(r'$z$')
-    ax.set_ylabel(r'$m_r$')
+    ax.set_xlabel(r'Redshift ($z$)')
+    ax.set_ylabel(r'Apparent Magnitude ($m_r$)')
     ax.set_xlim(xmnx)
     ax.set_ylim(ymnx)
 
-    ax.xaxis.set_major_locator(plt.MultipleLocator(0.5))
+    if logz_scale == True:
+        ax.set_xscale('log')
+        ax.set_xlim(1e-2,z_max+0.5)
+        ax.xaxis.set_major_locator(plt.LogLocator(base=10, numticks=12))
+    else:
+        ax.xaxis.set_major_locator(plt.MultipleLocator(0.5))
     plt.colorbar(c,label='p(z|DM) [a.u.]',ax=ax)
 
-    ax.legend()
-
+    ax.legend(loc='upper right', frameon=True, shadow=True)
+    # set the title of the figure
+    ax.set_title(title)
     frb_fig_u.set_fontsize(ax, 15.)
 
     # End

frb/galaxies/utils.py

@@ -6,6 +6,7 @@
 from pkg_resources import resource_filename
 import numpy as np
 from scipy.interpolate import interp1d
+import warnings
 
 import pandas
 
@@ -20,7 +21,7 @@
 from astropy import units
 
 import pandas as pd
-import extinction
+
 from linetools.spectra import xspectrum1d
 
 from frb import frb
@@ -37,6 +38,10 @@ def deredden_spec(spectrum:xspectrum1d.XSpectrum1D, ebv:float):
     """
 
     # Correct for Galactic extinction
+    # Hidnig this here to avoid a hard dependency
+    # TODO
+    #   Need to replace it 
+    import extinction
     AV = ebv * 3.1  # RV
     Al = extinction.ccm89(spectrum.wavelength.value, AV, 3.1)
     # New spec
@@ -70,7 +75,7 @@ def load_specdb(specdb_file=None):
             specdb_file = specdb_files[0]
             print("Loading spectra from {:s}".format(specdb_file))
         else:
-            raise IOError("There are no FRB_specdb.hdf5 files in your SPECDB folder")
+            raise IOError("There are no FRB_specDB_*.hdf5 files in your SPECDB folder")
     # Load it up
     specDB = SpecDB(db_file=specdb_file)
     # Return