orca.utils.beam_nmahesh ======================= .. py:module:: orca.utils.beam_nmahesh .. autoapi-nested-parse:: Primary beam models for OVRO-LWA (N. Mahesh implementation). This module provides classes for computing primary beam corrections using various beam models including analytic approximations and simulated beam patterns. Used for source flux density corrections and calibrator modeling. Classes: analytic_beam: Simple sin^1.6 elevation-based beam model. jones_beam: Full Jones matrix beam from HDF5 simulation files. Classes ------- .. autoapisummary:: orca.utils.beam_nmahesh.analytic_beam orca.utils.beam_nmahesh.woody_beam orca.utils.beam_nmahesh.jones_beam Functions --------- .. autoapisummary:: orca.utils.beam_nmahesh.knn_search orca.utils.beam_nmahesh.primary_beam_correction_val Module Contents --------------- .. py:function:: knn_search(arr, grid) Find nearest neighbor of array of points in a multi-dimensional grid. :param arr: Query points array. :param grid: Grid of reference points. :returns: Index of nearest neighbor. Source: glowingpython.blogspot.com/2012/04/k-nearest-neighbor-search.html .. py:function:: primary_beam_correction_val(pol, jones_matrix) Calculate primary beam correction value for a polarization. :param pol: Polarization type ('XX', 'YY', or 'I'). :param jones_matrix: 2x2 Jones matrix. :returns: Beam correction factor for the specified polarization. .. py:class:: analytic_beam(msfile=None, beam_file_path='/opt/beam/', freq=None) Analytic primary beam model using sin^1.6 elevation dependence. A simplified beam model that approximates the OVRO-LWA primary beam as a function of elevation only, ignoring azimuthal structure. :param msfile: Path to measurement set (unused). :param beam_file_path: Path to beam files (unused). :param freq: Observing frequency (unused). .. py:method:: srcjones(az, el) .. py:method:: read_beam_file() dummy function to make the code similar for all classes .. py:method:: get_muller_matrix_XY(jones_matrix) The muller matrix is in XY coordinate frame So (XX_obs,XY_obs,YX_obs,YY_obs).T=((M00,M01,M02,M03), (M10,M11, M12, M13), (M20, M21, M22, M23), (M30, M31, M32, M33)) (XX,XY,YX,YY).T Note that this Muller Matrix is not normalised. .. py:method:: get_source_pol_factors(jones_matrix) :staticmethod: .. py:method:: get_muller_matrix_stokes(jones_matrix) The muller matrix is in Stokes coordinate frame So (I_obs,Q_obs,U_obs,V_obs).T=((M00,M01,M02,M03), (M10,M11, M12, M13), (M20, M21, M22, M23), (M30, M31, M32, M33)) (I,Q,U,V).T Note that this Muller Matrix is not normalised. .. py:class:: woody_beam(msfile=None, beam_file_path='/opt/beam/', freq=None) This class uses a beam_file_path. This path is the location where all the relevant beamfiles are located. If freq is provided, then msfile is not required. However, if freq is not provided, msfile is must. In this case, the frequency used is the central frequency corresponding to the msfile. This code can only accept one frequency, but can work simulatenously over multiple sources. This code uses the Woody beam model. .. py:attribute:: beamIQUV :value: None .. py:attribute:: Ibeam :value: None .. py:attribute:: Ubeam :value: None .. py:attribute:: Qbeam :value: None .. py:attribute:: Vbeam :value: None .. py:attribute:: freq :value: None .. py:method:: ctrl_freq() Reads the central frequency of the MS file for the spw 0. .. py:method:: get_beam_file() Gets all the beamfiles from the relevant location and then chooses the file while is most relevant for the frequency which the user has requested. .. py:method:: read_beam_file() Reads the beamfile. If reading the beamfile is unsuccessfull, then it switches to the analytical beam. .. py:property:: beam_file_path .. py:property:: msfile .. py:method:: srcjones(az, el) The function name is srcjones to keep naming consistent with Jones_beam class. But this function returns I,Q,U,V beams. Compute beam scaling factor Can work simulatenously over multiple sources. So the az and el can be two arrays. The number of sources is assumed to be equal to the number of elements in alt/az :param az: azimuth in degrees :param el: elevation in degrees Returns: [I,Q,U,V] flux factors, where for an unpolarized source [I,Q,U,V] = [1,0,0,0] .. py:method:: get_source_pol_factors(jones_matrix) :staticmethod: I am assuming that the source is unpolarised. At these low frequencies this is a good assumption. Since the jones matrix in this class is essentially the source pol factors, just returning. .. py:method:: get_muller_matrix_XY(jones_matrix) The muller matrix is in XY coordinate frame So (XX_obs,XY_obs,YX_obs,YY_obs).T=((M00,M01,M02,M03), (M10,M11, M12, M13), (M20, M21, M22, M23), (M30, M31, M32, M33)) (XX,XY,YX,YY).T Note that this Muller Matrix is not normalised. .. py:method:: get_muller_matrix_stokes(jones_matrix) The muller matrix is in Stokes coordinate frame So (I_obs,Q_obs,U_obs,V_obs).T=((M00,M01,M02,M03), (M10,M11, M12, M13), (M20, M21, M22, M23), (M30, M31, M32, M33)) (I,Q,U,V).T Note that this Muller Matrix is not normalised. .. py:class:: jones_beam(beam_file_path='/lustre/msurajit/beam_model_nivedita/OVRO-LWA_soil_pt.h5', msfile=None, freq=None) For loading and returning LWA dipole beam values (derived from simulations made by Nivedita) Can only take one frequency at a time now. The beamfile path should contain the absolute path of the beamfile. If freq is provided, then msfile is not required. However, if freq is not provided, msfile is must. In this case, the frequency used is the central frequency corresponding to the msfile. .. py:property:: beamfile .. py:attribute:: freq :value: None .. py:attribute:: num_freq .. py:method:: ctrl_freq() Reads the central frequency of the MS file for the spw 0. .. py:property:: msfile .. py:method:: read_beam_file() az,za units: radian. Note that while the beamfiles accept in radian units, this function takes the alt,azimuth in degrees. The altitude is converted to zenith angle within the code. The jones matrices will be normalised with respect to the zenith .. py:method:: match_dimensions(data) .. py:method:: get_max_val() does accurate beam normalisation .. py:method:: srcjones(az, el) Compute beam scaling factor Can work simultaneously over multiple sources. So the az and el can be two arrays. The number of sources is assumed to be equal to the number of elements in alt/az We have implemented an approximate normalization of the beam, which is used by default. While the function to do a more accurate normalization is already implemented, it is not used by default. :param (az: :param el) coordinates in degrees: Returns: Jones matrix at coordinates (az,el) .. py:method:: get_source_pol_factors(jones_matrix) I am assuming that the source is unpolarised. At these low frequencies this is a good assumption. .. py:method:: get_muller_matrix_XY(jones_matrix) The muller matrix is in XY coordinate frame So (XX_obs,XY_obs,YX_obs,YY_obs).T=((M00,M01,M02,M03), (M10,M11, M12, M13), (M20, M21, M22, M23), (M30, M31, M32, M33)) (XX,XY,YX,YY).T Note that this Muller Matrix is not normalised. .. py:method:: get_muller_matrix_stokes(jones_matrix) The muller matrix is in Stokes coordinate frame So (I_obs,Q_obs,U_obs,V_obs).T=((M00,M01,M02,M03), (M10,M11, M12, M13), (M20, M21, M22, M23), (M30, M31, M32, M33)) (I,Q,U,V).T Note that this Muller Matrix is not normalised.