jit added

datamodel.rst

@@ -0,0 +1,62 @@
+qsoabsfind
+==========
+
+Instructions
+-------------
+
+Please read and perform the following steps.
+
+Running Help Command
+--------------------
+
+To get an overview of the tool's functionality and options, run:
+
+::
+
+    qsoabsfind --help
+
+Input FITS File Structure
+-------------------------
+
+The input `fits file` must have the following HDU extensions:
+
+- **FLUX**: Should ideally contain the residual spectra (usually the flux/continuum, i.e., the continuum normalized spectra).
+- **WAVELENGTH**: Observed wavelength (in Angstroms).
+- **ERROR**: Error on residuals.
+- **TGTDETALS**: Spectral details (such as Z_QSO, RA_QSO, DEC_QSO).
+
+Constant File (Optional)
+------------------------
+
+The user-defined **constant-file** must follow the same structure as the `qsoabsfind.constants` file, otherwise, the code will fail. If you want to use the default search parameters, you can run the tool without the `constant-file` option.
+
+Running the Tool
+----------------
+
+Run `qsoabsfind` with the required FITS file. If using a custom constant file, include it in the command:
+
+::
+
+    qsoabsfind --input <path_to_input_fits_file> [--constant-file <path_to_constant_file>] --output <path_to_output_fits_file>
+
+Replace `<path_to_input_fits_file>` with the path to your input FITS file, `<path_to_constant_file>` with the path to your constant file (if using), and `<path_to_output_fits_file>` with the desired path for your output FITS file.
+
+Output FITS File Structure
+--------------------------
+
+The **output** `fits file` will have the `ABSORBER` HDU, containing arrays such as:
+
+- **INDEX_SPEC**: Index of quasar (can be used to read the RA, DEC, and Z of QSOs).
+- **Z_ABS**: Redshift of absorber.
+- **${METAL}_${line}_EW**: Rest-frame equivalent widths (EWs) of absorber lines (e.g., MgII 2796, 2803 or CIV 1548, 1550) in Angstroms.
+- **${METAL}_${line}_EW_ERROR**: Uncertainties in rest-frame EWs of absorber lines in Angstroms.
+- **Z_ABS_ERR**: Measured error in the redshift of the absorber.
+- **GAUSS_FIT**: Rest-frame fitting parameters of double Gaussian to the absorber doublet (the width can be used to measure the velocity dispersion).
+- **GAUSS_FIT_STD**: Uncertainties in rest-frame fitting parameters of double Gaussian to the absorber doublet.
+- **SN_${METAL}_${line}**: Signal-to-noise ratio of the lines.
+- **${metal}_EW_TOTAL**: Total EW of the lines in Angstroms.
+- **${metal}_EW_TOTAL_ERROR**: Uncertainties in total EW of the lines in Angstroms.
+
+Thanks,
+Abhijeet Anand
+Lawrence Berkeley National Lab  

qsoabsfind/absfinder.py

@@ -1,7 +1,7 @@
 import numpy as np
 from functools import reduce
 from operator import add
-from .utils import convolution_fun, save_plot, plot_absorber
+from .utils import convolution_fun
 from .absorberutils import (
     estimate_local_sigma_conv_array, group_and_weighted_mean_selection_function,
     median_selection_after_combining,
@@ -11,7 +11,9 @@
 from .ew import measure_absorber_properties_double_gaussian
 from .config import lines, speed_of_light, oscillator_parameters
 from .spec import QSOSpecRead
+from numba import jit
 
+@jit(nopython=True)
 def find_valid_indices(our_z, residual_our_z, lam_search, conv_arr, sigma_cr, coeff_sigma, d_pix, beta, line1, line2):
     """
     Find valid indices based on thresholding in the convolved array.
@@ -139,17 +141,14 @@ def convolution_method_absorber_finder_in_QSO_spectra(fits_file, spec_index, abs
             sigma_cr = estimate_local_sigma_conv_array(conv_arr, pm_pixel=pm_pixel)
             thr = np.nanmedian(conv_arr) - coeff_sigma * sigma_cr
 
-            if spec_index%20==0:
-                _ = convolution_fun(absorber, mult_resi * unmsk_residual, sig_ker, amp_ratio=0.5, log=logwave, index=spec_index)
-                save_plot(lam_search, [unmsk_residual, conv_arr, thr], f'qsoabsfind_plot_spec_{spec_index}.png', 'wave', 'flux', f'{spec_index}')
 
             conv_arr[np.isnan(conv_arr)] = 1e5
             our_z_ind = conv_arr < thr
             conv_arr[conv_arr == 1e5] = np.nan
 
             our_z = lam_search[our_z_ind] / line_centre - 1
             residual_our_z = unmsk_residual[our_z_ind]
-            #print(f'INFO:: Potential absorber = {our_z}')
+
             new_our_z, new_res_arr = find_valid_indices(our_z, residual_our_z, lam_search, conv_arr, sigma_cr, coeff_sigma, d_pix, f1 / f2, line1, line2)
 
             final_our_z = group_and_weighted_mean_selection_function(new_our_z, np.array(new_res_arr))
@@ -165,15 +164,11 @@ def convolution_method_absorber_finder_in_QSO_spectra(fits_file, spec_index, abs
         combined_final_our_z = combined_final_our_z[~np.isnan(combined_final_our_z)]
         combined_final_our_z = combined_final_our_z.tolist()
 
-        #plot_absorber(lam_obs, residual, absorber, combined_final_our_z, xlabel='obs wave (ang)', ylabel='residual', title=f'QSO_{spec_index}', plot_filename=None)
-
         if len(combined_final_our_z)>0:
 
             z_abs, z_err, fit_param, fit_param_std, EW_first_line_mean, EW_second_line_mean, EW_total_mean, EW_first_line_error, EW_second_line_error, EW_total_error = measure_absorber_properties_double_gaussian(
                 index=spec_index, wavelength=lam_obs, flux=residual, error=error, absorber_redshift=combined_final_our_z, bound=bound, use_kernel=absorber, d_pix=d_pix, use_covariance=use_covariance)
 
-
-            #print(f'INFO:: {spec_index}, {z_abs}, {fit_param}')
             pure_z_abs = np.zeros(len(z_abs))
             pure_gauss_fit = np.zeros((len(z_abs), 6))
             pure_gauss_fit_std = np.zeros((len(z_abs), 6))
@@ -240,7 +235,7 @@ def convolution_method_absorber_finder_in_QSO_spectra(fits_file, spec_index, abs
                 if absorber=='MgII':
                     match_abs1 = remove_Mg_falsely_identified_as_Fe_absorber(spec_index, pure_z_abs, lam_obs, residual, error, d_pix=d_pix)
                 else:
-                    match_abs1 = -1
+                    match_abs1 = -1*np.ones(len(pure_z_abs))
                 match_abs2 = z_abs_from_same_metal_absorber(pure_z_abs, lam_obs, residual, error, d_pix, use_kernel=absorber)
                 ind_z = contiguous_pixel_remover(pure_z_abs, sn1_all, sn2_all)
                 sel_indices = (match_abs1 == -1) & (match_abs2 == -1) & (ind_z == -1)  # pure final absorber candidates

qsoabsfind/absorberutils.py

@@ -9,8 +9,7 @@
 from .utils import elapsed
 from numba import jit
 
-
-
+@jit(nopython=True)
 def estimate_local_sigma_conv_array(conv_array, pm_pixel):
     """
     Estimate the local standard deviation for each element in the convolution array over a window defined by pm_pixel.
@@ -380,6 +379,7 @@ def contiguous_pixel_remover(abs_z, sn1_all, sn2_all, use_kernel='MgII'):
 
     return ind_true
 
+@jit(nopython=True)
 def redshift_estimate(fitted_obs_l1, fitted_obs_l2, std_fitted_obs_l1, std_fitted_obs_l2, line1, line2):
     """
     Estimate the redshift and its error from Gaussian fitting parameters for two spectral lines.

qsoabsfind/utils.py

@@ -87,9 +87,6 @@ def convolution_fun(absorber, residual_arr_after_mask, width, amp_ratio=0.5, log
 
     gauss_kernel = gauss_two_lines_kernel(lam_ker, a=ker_parm)
 
-    if index is not None:
-        save_plot(lam_ker, gauss_kernel, f'kernel_{index}.png', 'wave', 'flux', f'{index}')
-
     result = np.convolve(gauss_kernel, residual_arr_after_mask, mode='same')
 
     #check if input and output array size are same