Source code for irfpy.cena.fov

r''' fov module is for Chap 3. in the CENA FM calibration report.

The :ref:`CENA calibration report issue 1 rev 1 <cenacalrep>` is implemented here.
The module treats Chapter 3: Sensor properties.

A simple use of the module can be found as a script of :mod:`cena_visualize_fov`.

Conversion between |theta| |phi| system to a vector is implemented in :mod:`irfpy.cena.frame` module.
See also :ref:`frame_conversion`.

Definition of |theta| and |phi| is as follows.
Primary axis is -*z* axis.
|theta| is measured from *xy* plane and positive in -*z* axis.
In other word, this is a co-latitude measured from -*z* axis.
|phi| is measured from *y* axis (strictly speaking *yz* plane) and
positive toward *x* axis.
See Figure 1.1 in :ref:`cenacalrep`.
'''
import math
import numpy as np


def shape_of_response(n, theta, phi):
    ''' Return the shape of response.

    |theta| is elevation angle and |phi| is the azimuthal angle in degrees.
    |theta| should be -90 to 90, and |phi| should be -180 to 180.
    *n* is the channel ID from 0 to 6.
    '''
    if n < 0 or n > 6:
        raise IndexError('*n* should be between 0 and 6')
    if theta < -90 or theta > 90:
        raise IndexError('*theta* should be between -90 and 90')
    if phi < -180 or phi > 180:
        raise IndexError('*phi* should be between -180 and 180')

    sig = azim_pix_fwhm(n) / 2.35
    yet = elev_pix_fwhm(n) / 2.35

    phi_c = azim_pix_center(n)
    theta_c = elev_pix_center(n)

    coeff = - (((phi - phi_c) ** 2) / (2 * sig * sig)
                + ((theta - theta_c) ** 2) / (2 * yet * yet))
    return math.exp(coeff)

def pixel_shape(n, elevresolution=0, azimresolution=0, resolution=None, fill=False):
    r''' Return an array of the pixel edge angle of its FWHM.

    :param n: Channel number.  0 to 6.
    :param azimresolution: Azimuthal resolution. Number of division in azimuthal direction.
    :param elevresolution: Elevation resolution. Number of division in elevation direction.
    :param resolution: Resolution. Overwrite ``azimresolution`` and ``elevresolution``.
    :param fill: If true, the first point is added to the end.
    :returns: A tuple of 2 numpy array. Array of :math:`\theta` and array of :math:`\phi`.

    The number of elements is ``4+2*azimresolution+2*elevresolution``.
    If fill is True, one more element is added.
    However, it is recommend to check the number after you get the array.

    >>> pix0 = pixel_shape(0, azimresolution=50, elevresolution=3)
    >>> print(len(pix0[0]), len(pix0[1]))
    110 110

    '''
    if resolution is not None:
        azimresolution = resolution
        elevresolution = resolution
    azimc = azim_pix_center(n)
    elevc = elev_pix_center(n)
    azimw = azim_pix_fwhm(n)
    elevw = elev_pix_fwhm(n)

    azim0 = azimc - azimw / 2.
    elev0 = elevc - elevw / 2.

    azim1 = azimc + azimw / 2.
    elev1 = elevc + elevw / 2.

    azndiv = azimresolution + 1
    dazim = azimw / azndiv
    elndiv = elevresolution + 1
    delev = elevw / elndiv

    az_list = []
    el_list = []

    for idiv in range(elndiv):
        az_list.append(azim1)
        el_list.append(elev0 + idiv * delev)
    for idiv in range(azndiv):
        az_list.append(azim1 - idiv * dazim)
        el_list.append(elev1)
    for idiv in range(elndiv):
        az_list.append(azim0)
        el_list.append(elev1 - idiv * delev)
    for idiv in range(azndiv):
        az_list.append(azim0 + idiv * dazim)
        el_list.append(elev0)

    if fill:
        az_list.append(az_list[0])
        el_list.append(el_list[0])

    return np.array(el_list), np.array(az_list)
        

def azim_pix_center(n):
    ''' Returns the azimuth (|phi|) center in degrees for the pixel *n*.

    Implements equation (3.2), -19.06*n + 57.19
    '''
    if n < 0 or n > 6:
        raise IndexError('*n* should be between 0 and 6')
    return -19.06 * n + 57.19


def azim_pix_fwhm(n):
    ''' Returns the FWHM of the azimuthal directions (|phi|) in degrees.

    Implements equation (3.3).  Indeed it does not depend on *n*,
    45.0 degrees.
    '''
    return 45.0


def elev_pix_center(n):
    ''' Returns the elevation center (|theta|) in degrees for the pixel *n*.

    Implements equation (3.4). It does not indeed depend on *n*:
    -6.05 degrees.
    '''
    return -6.05


def elev_pix_fwhm(n):
    ''' Returns the FWHM of the elevation (|theta|) directions in degrees.

    Implements equation (3.5).  Independent of *n*, 6.44 is returned.
    '''
    return 6.44


def solid_angle(n):
    r''' Returns the approximate soliad angle.

    Equation (3.8) gives the definition.

    .. math::

        \Delta\Omega_n \approx FWHM_{\phi,n} \cdot FWHM_{\theta,n}

    >>> print('%.3f' % solid_angle(0))
    0.088
    '''
    elev = elev_pix_fwhm(n) * math.pi / 180.
    azim = azim_pix_fwhm(n) * math.pi / 180.
    return elev * azim


def infov(theta, phi):
    r""" A simple FOV envelope checker.

    :param theta: The |theta| angle(s), from -90 to 90 in degrees. Looking direction.
    :param phi: The |phi| angle(s), from -180 to 180 in degrees. Looking direction.

    >>> infov(0, 0)
    True
    >>> infov(5, 0)    # Theta=5, False
    False
    >>> infov(-5, 0)
    True
    >>> infov([-6, -6, -6, -6, -6], [-180, -90, 0, 90, 180])
    array([False, False,  True, False, False])


    """
    _theta = np.array(theta)
    _phi = np.array(phi)
    theta_cond = np.logical_and(_theta >= -15, _theta <= 0)
    phi_cond = np.logical_and(_phi >= -75, _phi <= 75)

    return np.logical_and(theta_cond, phi_cond)