rm_lite.utils.synthesis ======================= .. py:module:: rm_lite.utils.synthesis .. autoapi-nested-parse:: RM-synthesis utils Attributes ---------- .. autoapisummary:: rm_lite.utils.synthesis.T rm_lite.utils.synthesis.fdf_params_schema rm_lite.utils.synthesis.fdf_params_schema_df Classes ------- .. autoapisummary:: rm_lite.utils.synthesis.FDFOptions rm_lite.utils.synthesis.FWHM rm_lite.utils.synthesis.FractionalSpectra rm_lite.utils.synthesis.RMSFParams rm_lite.utils.synthesis.RMSFResults rm_lite.utils.synthesis.RMSynthParams rm_lite.utils.synthesis.RMsynthResults rm_lite.utils.synthesis.SigmaAdd rm_lite.utils.synthesis.SigmaAddArrays rm_lite.utils.synthesis.StokesData rm_lite.utils.synthesis.StokesSigmaAdd rm_lite.utils.synthesis.TheoreticalNoise Functions --------- .. autoapisummary:: rm_lite.utils.synthesis.calc_mom2_fdf rm_lite.utils.synthesis.calculate_sigma_add rm_lite.utils.synthesis.calculate_sigma_add_arr rm_lite.utils.synthesis.cdf_percentile rm_lite.utils.synthesis.compute_rmsf_params rm_lite.utils.synthesis.compute_rmsynth_params rm_lite.utils.synthesis.compute_theoretical_noise rm_lite.utils.synthesis.create_fractional_spectra rm_lite.utils.synthesis.faraday_simple_spectrum rm_lite.utils.synthesis.freq_to_lambda2 rm_lite.utils.synthesis.get_fdf_parameters rm_lite.utils.synthesis.get_fwhm_rmsf rm_lite.utils.synthesis.get_mask_index rm_lite.utils.synthesis.get_rmsf_nufft rm_lite.utils.synthesis.inverse_rmsynth_nufft rm_lite.utils.synthesis.lambda2_to_freq rm_lite.utils.synthesis.make_double_phi_arr rm_lite.utils.synthesis.make_phi_arr rm_lite.utils.synthesis.measure_fdf_complexity rm_lite.utils.synthesis.measure_qu_complexity rm_lite.utils.synthesis.rmsynth_nufft Module Contents --------------- .. py:class:: FDFOptions Bases: :py:obj:`NamedTuple` Options for RM-synthesis .. py:attribute:: d_phi_radm2 :type: float | None :value: None Faraday depth resolution .. py:attribute:: do_fit_rmsf :type: bool :value: False Fit RMSF .. py:attribute:: do_fit_rmsf_real :type: bool :value: False Fit real part of the RMSF .. py:attribute:: n_samples :type: float | None :value: 10.0 Number of samples .. py:attribute:: phi_max_radm2 :type: float | None :value: None Maximum Faraday depth .. py:attribute:: weight_type :type: Literal['variance', 'uniform'] :value: 'variance' Weight type .. py:class:: FWHM Bases: :py:obj:`NamedTuple` .. py:attribute:: d_lambda_sq_max_m2 :type: float The maximum difference in lambda^2 values .. py:attribute:: fwhm_rmsf_radm2 :type: float The FWHM of the RMSF main lobe .. py:attribute:: lambda_sq_range_m2 :type: float The range of lambda^2 values .. py:class:: FractionalSpectra Bases: :py:obj:`NamedTuple` .. py:attribute:: fit_result :type: rm_lite.utils.fitting.FitResult | None .. py:attribute:: no_nan_idx :type: numpy.typing.NDArray[numpy.bool_] .. py:attribute:: stokes_data :type: StokesData .. py:class:: RMSFParams Bases: :py:obj:`NamedTuple` RM spread function parameters .. py:attribute:: phi_max :type: float Maximum Faraday depth .. py:attribute:: phi_max_scale :type: float Maximum Faraday depth scale .. py:attribute:: rmsf_fwhm_meas :type: float Measured FWHM of the RMSF .. py:attribute:: rmsf_fwhm_theory :type: float Theoretical FWHM of the RMSF .. py:class:: RMSFResults Bases: :py:obj:`NamedTuple` Results of the RMSF calculation .. py:attribute:: fit_status_arr :type: numpy.typing.NDArray[numpy.float64] The status of the RMSF fit .. py:attribute:: fwhm_rmsf_arr :type: numpy.typing.NDArray[numpy.float64] The FWHM of the RMSF main lobe .. py:attribute:: phi_double_arr_radm2 :type: numpy.typing.NDArray[numpy.float64] The (double length) Faraday depth array .. py:attribute:: rmsf_cube :type: numpy.typing.NDArray[numpy.float64] The RMSF cube .. py:class:: RMSynthParams Bases: :py:obj:`NamedTuple` Parameters for RM-synthesis calculation .. py:attribute:: lam_sq_0_m2 :type: float Reference wavelength^2 value .. py:attribute:: lambda_sq_arr_m2 :type: numpy.typing.NDArray[numpy.float64] Wavelength^2 values in m^2 .. py:attribute:: phi_arr_radm2 :type: numpy.typing.NDArray[numpy.float64] Faraday depth values in rad/m^2 .. py:attribute:: weight_arr :type: numpy.typing.NDArray[numpy.float64] Weight array .. py:class:: RMsynthResults Bases: :py:obj:`NamedTuple` Results of the RM-synthesis calculation .. py:attribute:: fdf_dirty_cube :type: numpy.typing.NDArray[numpy.float64] The Faraday dispersion function cube .. py:attribute:: lam_sq_0_m2 :type: float The reference lambda^2 value .. py:class:: SigmaAdd Bases: :py:obj:`NamedTuple` Sigma_add complexity metrics .. py:attribute:: sigma_add :type: float Sigma_add median value .. py:attribute:: sigma_add_arrays :type: SigmaAddArrays Sigma_add arrays .. py:attribute:: sigma_add_minus :type: float Sigma_add lower quartile .. py:attribute:: sigma_add_plus :type: float Sigma_add upper quartile .. py:class:: SigmaAddArrays Bases: :py:obj:`NamedTuple` .. py:attribute:: cdf :type: numpy.typing.NDArray[numpy.float64] CDF array of the additional noise term .. py:attribute:: pdf :type: numpy.typing.NDArray[numpy.float64] PDF array of the additional noise term .. py:attribute:: sigma_add_arr :type: numpy.typing.NDArray[numpy.float64] Array of additional noise values .. py:class:: StokesData Bases: :py:obj:`NamedTuple` Stokes parameters and errors .. py:method:: with_options(**kwargs) Create a new StokesData instance with keywords updated :returns: New StokesData instance with updated attributes :rtype: StokesData .. py:attribute:: complex_pol_arr :type: numpy.typing.NDArray[numpy.complex128] Stokes Q and U array .. py:attribute:: complex_pol_error :type: numpy.typing.NDArray[numpy.complex128] Stokes Q and U error array .. py:attribute:: freq_arr_hz :type: numpy.typing.NDArray[numpy.float64] Frequency array in Hz .. py:attribute:: stokes_i_arr :type: numpy.typing.NDArray[numpy.float64] | None :value: None Stokes I array .. py:attribute:: stokes_i_error_arr :type: numpy.typing.NDArray[numpy.float64] | None :value: None Stokes I error array .. py:attribute:: stokes_i_model_arr :type: numpy.typing.NDArray[numpy.float64] | None :value: None Stokes I model array .. py:attribute:: stokes_i_model_error :type: numpy.typing.NDArray[numpy.float64] | None :value: None Stokes I model error array .. py:class:: StokesSigmaAdd Bases: :py:obj:`NamedTuple` Stokes Sigma_add complexity metrics .. py:attribute:: sigma_add_p :type: SigmaAdd Sigma_add for polarised intensity .. py:attribute:: sigma_add_q :type: SigmaAdd Sigma_add for Stokes Q .. py:attribute:: sigma_add_u :type: SigmaAdd Sigma_add for Stokes U .. py:class:: TheoreticalNoise Bases: :py:obj:`NamedTuple` Theoretical noise of the FDF .. py:method:: with_options(**kwargs) Create a new TheoreticalNoise instance with keywords updated :returns: New TheoreticalNoise instance with updated attributes :rtype: TheoreticalNoise .. py:attribute:: fdf_error_noise :type: float Theoretical noise of the FDF .. py:attribute:: fdf_q_noise :type: float Theoretical noise of the real FDF .. py:attribute:: fdf_u_noise :type: float Theoretical noise of the imaginary FDF .. py:function:: calc_mom2_fdf(complex_fdf_arr: numpy.typing.NDArray[numpy.complex128], phi_arr_radm2: numpy.typing.NDArray[numpy.float64]) -> float Calculate the 2nd moment of the polarised intensity FDF. Can be applied to a clean component spectrum or a standard FDF .. py:function:: calculate_sigma_add(y_arr: numpy.typing.NDArray[numpy.float64], dy_arr: numpy.typing.NDArray[numpy.float64], median: float | None = None, noise: float | None = None, n_samples: int = 1000) -> SigmaAdd Calculate the most likely value of additional scatter, assuming the input data is drawn from a normal distribution. The total uncertainty on each data point Y_i is modelled as dYtot_i**2 = dY_i**2 + dYadd**2. .. py:function:: calculate_sigma_add_arr(y_arr: numpy.typing.NDArray[numpy.float64], dy_arr: numpy.typing.NDArray[numpy.float64], median: float | None = None, noise: float | None = None, n_samples: int = 1000) -> SigmaAddArrays .. py:function:: cdf_percentile(values: numpy.typing.NDArray[numpy.float64], cdf: numpy.typing.NDArray[numpy.float64], q=50.0) -> float Return the value at a given percentile of a cumulative distribution function :param values: Array of values :type values: NDArray[np.float64] :param cdf: Cumulative distribution function :type cdf: NDArray[np.float64] :param q: Percentile. Defaults to 50.0. :type q: float, optional :returns: Interpolated value at the given percentile :rtype: float .. py:function:: compute_rmsf_params(freq_arr_hz: numpy.typing.NDArray[numpy.float64], weight_arr: numpy.typing.NDArray[numpy.float64]) -> RMSFParams .. py:function:: compute_rmsynth_params(freq_arr_hz: numpy.typing.NDArray[numpy.float64], complex_pol_arr: numpy.typing.NDArray[numpy.complex128], complex_pol_error: numpy.typing.NDArray[numpy.complex128], fdf_options: FDFOptions) -> RMSynthParams Calculate the parameters for RM-synthesis. :param freq_arr_hz: Frequency array in Hz :type freq_arr_hz: NDArray[np.float64] :param pol_arr: Complex polarisation array :type pol_arr: NDArray[np.complex128] :param real_qu_error: Error in Stokes Q and U (real) :type real_qu_error: NDArray[np.float64 | np.float32] :param fdf_options: Options for RM-synthesis :type fdf_options: FDFOptions :raises ValueError: If d_phi_radm2 is not provided and n_samples is None. :returns: Wavelength^2 values, reference wavelength^2, Faraday depth values, weight array :rtype: RMSynthParams .. py:function:: compute_theoretical_noise(complex_pol_error: numpy.typing.NDArray[numpy.complex128], weight_arr: numpy.typing.NDArray[numpy.float64]) -> TheoreticalNoise .. py:function:: create_fractional_spectra(stokes_data: StokesData, ref_freq_hz: float, fit_order: int = 2, fit_function: Literal['log', 'linear'] = 'log', n_error_samples: int = 10000) -> FractionalSpectra | None .. py:function:: faraday_simple_spectrum(freq_arr_hz: numpy.typing.NDArray[numpy.float64], frac_pol: float, psi0_deg: float, rm_radm2: float) -> numpy.typing.NDArray[numpy.complex128] Create a simple Faraday spectrum with a single component. :param freq_arr_hz: Frequency array in Hz :type freq_arr_hz: NDArray[np.float64] :param frac_pol: Fractional polarization :type frac_pol: float :param psi0_deg: Initial polarization angle in degrees :type psi0_deg: float :param rm_radm2: RM in rad/m^2 :type rm_radm2: float :returns: Complex polarization spectrum :rtype: NDArray[np.float64] .. py:function:: freq_to_lambda2(freq_hz: T) -> T Convert frequency to lambda^2. :param freq_hz: Frequency in Hz :type freq_hz: float :returns: Wavelength^2 in m^2 :rtype: float .. py:function:: get_fdf_parameters(fdf_arr: numpy.typing.NDArray[numpy.complex128], phi_arr_radm2: numpy.typing.NDArray[numpy.float64], fwhm_rmsf_radm2: float, freq_arr_hz: numpy.typing.NDArray[numpy.float64], complex_pol_arr: numpy.typing.NDArray[numpy.complex128], complex_pol_error: numpy.typing.NDArray[numpy.complex128], lambda_sq_arr_m2: numpy.typing.NDArray[numpy.float64], lam_sq_0_m2: float, stokes_i_reference_flux: float, theoretical_noise: TheoreticalNoise, fit_function: Literal['log', 'linear'], bias_correction_snr: float = 5.0) -> polars.DataFrame Measure standard parameters from a complex Faraday Dispersion Function. Currently this function assumes that the noise levels in the Stokes Q and U spectra are the same. Returns a dictionary containing measured parameters. .. py:function:: get_fwhm_rmsf(lambda_sq_arr_m2: numpy.typing.NDArray[numpy.float64]) -> FWHM Calculate the FWHM of the RMSF. :param lambda_sq_arr_m2: Wavelength^2 values in m^2 :type lambda_sq_arr_m2: NDArray[np.float64] :param super_resolution: Use Cotton+Rudnick superresolution. Defaults to False. :type super_resolution: bool, optional :returns: FWHM of the RMSF main lobe, maximum difference in lambda^2 values, range of lambda^2 values :rtype: fwhm_rmsf_arr .. py:function:: get_mask_index(stokes_data: StokesData) -> numpy.typing.NDArray[numpy.bool_] .. py:function:: get_rmsf_nufft(lambda_sq_arr_m2: numpy.typing.NDArray[numpy.float64], phi_arr_radm2: numpy.typing.NDArray[numpy.float64], weight_arr: numpy.typing.NDArray[numpy.float64], lam_sq_0_m2: float, mask_arr: numpy.typing.NDArray[numpy.bool_] | None = None, do_fit_rmsf: bool = False, do_fit_rmsf_real=False, eps: float = 1e-06) -> RMSFResults Compute the RMSF for a given set of lambda^2 values. :param lambda_sq_arr_m2: Wavelength^2 values in m^2 :type lambda_sq_arr_m2: NDArray[np.float64] :param phi_arr_radm2: Faraday depth values in rad/m^2 :type phi_arr_radm2: NDArray[np.float64] :param weight_arr: Weight array :type weight_arr: NDArray[np.float64] :param lam_sq_0_m2: Reference wavelength^2 value :type lam_sq_0_m2: float :param super_resolution: Use superresolution. Defaults to False. :type super_resolution: bool, optional :param mask_arr: Mask array. Defaults to None. :type mask_arr: Optional[NDArray[np.float64]], optional :param do_fit_rmsf: Fit the RMSF with a Gaussian. Defaults to False. :type do_fit_rmsf: bool, optional :param do_fit_rmsf_real: Fit the *real* part of the. Defaults to False. :type do_fit_rmsf_real: bool, optional :param eps: NUFFT tolerance. Defaults to 1e-6. :type eps: float, optional :raises ValueError: If the wavelength^2 and weight arrays are not the same shape. :raises ValueError: If the mask dimensions are > 3. :raises ValueError: If the mask depth does not match the lambda^2 vector. :returns: rmsf_cube, phi_double_arr_radm2, fwhm_rmsf_arr, fit_status_arr :rtype: RMSFResults .. py:function:: inverse_rmsynth_nufft(complex_fdf_arr: numpy.typing.NDArray[numpy.complex128], lambda_sq_arr_m2: numpy.typing.NDArray[numpy.float64], phi_arr_radm2: numpy.typing.NDArray[numpy.float64], lam_sq_0_m2: float, eps: float = 1e-06) -> numpy.typing.NDArray[numpy.complex128] Inverse RM-synthesis - FDF to Stokes Q and U in wavelength^2 space. :param complex_fdf_arr: Complex polarisation array in Faraday depth space :type complex_fdf_arr: NDArray[np.complex128] :param lambda_sq_arr_m2: Wavelength^2 values in m^2 :type lambda_sq_arr_m2: NDArray[np.float64] :param phi_arr_radm2: Faraday depth values in rad/m^2 :type phi_arr_radm2: NDArray[np.float64] :param lam_sq_0_m2: Reference wavelength^2 value :type lam_sq_0_m2: float :param eps: NUFFT tolerance. Defaults to 1e-6. :type eps: float, optional :raises ValueError: If the Stokes Q and U data arrays are not the same shape. :raises ValueError: If the data dimensions are > 3. :raises ValueError: If the data depth does not match the lambda^2 vector. :returns: Complex polarisation array in wavelength^2 space :rtype: NDArray[np.float64] .. py:function:: lambda2_to_freq(lambda_sq_m2: T) -> T Convert lambda^2 to frequency. :param lambda_sq_m2: Wavelength^2 in m^2 :type lambda_sq_m2: NDArray[np.float64] :returns: Frequency in Hz :rtype: NDArray[np.float64] .. py:function:: make_double_phi_arr(phi_arr_radm2: numpy.typing.NDArray[numpy.float64]) -> numpy.typing.NDArray[numpy.float64] .. py:function:: make_phi_arr(phi_max_radm2: float, d_phi_radm2: float) -> numpy.typing.NDArray[numpy.float64] Construct a Faraday depth array. :param phi_max_radm2: Maximum Faraday depth in rad/m^2 :type phi_max_radm2: float :param d_phi_radm2: Spacing in Faraday depth in rad/m^2 :type d_phi_radm2: float :returns: Faraday depth array in rad/m^2 :rtype: NDArray[np.float64] .. py:function:: measure_fdf_complexity(phi_arr_radm2: numpy.typing.NDArray[numpy.float64], complex_fdf_arr: numpy.typing.NDArray[numpy.complex128]) -> float .. py:function:: measure_qu_complexity(freq_arr_hz: numpy.typing.NDArray[numpy.float64], complex_pol_arr: numpy.typing.NDArray[numpy.complex128], complex_pol_error: numpy.typing.NDArray[numpy.complex128], frac_pol: float, psi0_deg: float, rm_radm2: float) -> StokesSigmaAdd .. py:function:: rmsynth_nufft(complex_pol_arr: numpy.typing.NDArray[numpy.complex128], lambda_sq_arr_m2: numpy.typing.NDArray[numpy.float64], phi_arr_radm2: numpy.typing.NDArray[numpy.float64], weight_arr: numpy.typing.NDArray[numpy.float64], lam_sq_0_m2: float, eps: float = 1e-06) -> numpy.typing.NDArray[numpy.complex128] Run RM-synthesis on a cube of Stokes Q and U data using the NUFFT method. :param complex_pol_arr: Complex polarisation values (Q + iU) :type complex_pol_arr: NDArray[np.complex128] :param lambda_sq_arr_m2: Wavelength^2 values in m^2 :type lambda_sq_arr_m2: NDArray[np.float64] :param phi_arr_radm2: Faraday depth values in rad/m^2 :type phi_arr_radm2: NDArray[np.float64] :param weight_arr: Weight array :type weight_arr: NDArray[np.float64] :param lam_sq_0_m2: Reference wavelength^2 in m^2. Defaults to None. :type lam_sq_0_m2: Optional[float], optional :param eps: NUFFT tolerance. Defaults to 1e-6. :type eps: float, optional :raises ValueError: If the weight and lambda^2 arrays are not the same shape. :raises ValueError: If the Stokes Q and U data arrays are not the same shape. :raises ValueError: If the data dimensions are > 3. :raises ValueError: If the data depth does not match the lambda^2 vector. :returns: Dirty Faraday dispersion function cube :rtype: NDArray[np.float64] .. py:data:: T .. py:data:: fdf_params_schema .. py:data:: fdf_params_schema_df