Skip to content

Commit

Permalink
Add antenna_delays/weights and option to normalize or not the array_f…
Browse files Browse the repository at this point in the history
…actor
  • Loading branch information
AlanLoh committed Nov 27, 2023
1 parent ece2093 commit ef928ac
Show file tree
Hide file tree
Showing 3 changed files with 89 additions and 39 deletions.
2 changes: 1 addition & 1 deletion nenupy/__init__.py
Original file line number Diff line number Diff line change
Expand Up @@ -5,7 +5,7 @@
__copyright__ = "Copyright 2023, nenupy"
__credits__ = ["Alan Loh"]
__license__ = "MIT"
__version__ = "2.5.2"
__version__ = "2.5.3"
__maintainer__ = "Alan Loh"
__email__ = "alan.loh@obspm.fr"

Expand Down
92 changes: 68 additions & 24 deletions nenupy/instru/interferometer.py
Original file line number Diff line number Diff line change
Expand Up @@ -186,6 +186,8 @@ class Interferometer(ABC, metaclass=CombinedMeta):
~nenupy.instru.interferometer.Interferometer.baselines
~nenupy.instru.interferometer.Interferometer.position
~nenupy.instru.interferometer.Interferometer.size
~nenupy.instru.interferometer.Interferometer.antenna_weights
~nenupy.instru.interferometer.Interferometer.antenna_delays
.. rubric:: Methods Summary
Expand All @@ -209,15 +211,15 @@ def __init__(self,
antenna_names: np.ndarray,
antenna_positions: np.ndarray,
antenna_gains: np.ndarray,
extra_delays: np.ndarray = None,
feed_gains: np.ndarray = None,
antenna_delays: np.ndarray = None,
antenna_weights: np.ndarray = None,
):
self.position = position
self.antenna_names = antenna_names
self.antenna_positions = antenna_positions
self.antenna_gains = antenna_gains
self.extra_delays = extra_delays
self.feed_gains = feed_gains
self.antenna_delays = antenna_delays
self.antenna_weights = antenna_weights


def __getitem__(self, n):
Expand Down Expand Up @@ -259,6 +261,8 @@ def __getitem__(self, n):
interfero.antenna_names = self.antenna_names[antenna_mask]
interfero.antenna_positions = self.antenna_positions[antenna_mask, :]
interfero.antenna_gains = self.antenna_gains[antenna_mask]
interfero.antenna_weights = self.antenna_weights[antenna_mask]
interfero.antenna_delays = self.antenna_delays[antenna_mask]
return interfero


Expand All @@ -276,6 +280,8 @@ def __add__(self, other):
interferometer.antenna_names = np.concatenate((self.antenna_names, other.antenna_names[to_include]))
interferometer.antenna_positions = np.concatenate((self.antenna_positions, other.antenna_positions[to_include]))
interferometer.antenna_gains = np.concatenate((self.antenna_gains, other.antenna_gains[to_include]))
interferometer.antenna_weights = np.concatenate((self.antenna_weights, other.antenna_weights[to_include]))
interferometer.antenna_delays = np.concatenate((self.antenna_delays, other.antenna_delays[to_include]))
return interferometer


Expand All @@ -285,6 +291,8 @@ def __sub__(self, other):
interferometer.antenna_names = self.antenna_names[to_keep]
interferometer.antenna_positions = self.antenna_positions[to_keep]
interferometer.antenna_gains = self.antenna_gains[to_keep]
interferometer.antenna_weights = self.antenna_weights[to_keep]
interferometer.antenna_delays = self.antenna_delays[to_keep]
return interferometer


Expand Down Expand Up @@ -384,6 +392,42 @@ def position(self, p: EarthLocation):
self._position = p


@property
def antenna_weights(self) -> np.ndarray:
""" """
return self._antenna_weights
@antenna_weights.setter
def antenna_weights(self, weights: np.ndarray):
if weights is None:
weights = np.ones(self.size)
elif weights.size != self.size:
raise ValueError(
f"antenna_weights (of size {weights.size}) do not match antenna_positions (of shape {self.size})."
)
elif np.any(weights.min() < 0) or np.any(weights.max() > 1):
raise ValueError(
f"antenna_weights values are restricted between 0 and 1."
)
self._antenna_weights = weights


@property
def antenna_delays(self) -> np.ndarray:
""" Add delay errors for each antennae
They could be cable connection errors during construction, cables of wrong length, ...
"""
return self._antenna_delays
@antenna_delays.setter
def antenna_delays(self, delays: np.ndarray):
if delays is None:
delays = np.zeros(self.size)
elif delays.size != self.size:
raise ValueError(
f"antenna_delays (of size {delays.size}) do not match antenna_positions (of shape {self.size})."
)
self._antenna_delays = delays


# --------------------------------------------------------- #
# ------------------------ Methods ------------------------ #
def plot(self, **kwargs) -> None:
Expand Down Expand Up @@ -516,17 +560,17 @@ def plot(self, **kwargs) -> None:



def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = False) -> da.Array:
def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = False, normalize: bool = True) -> da.Array:
r""" Computes the array factor of the antenna distribution.
.. math::
\mathcal{F}(\nu, \phi, \theta) = \sum_{\rm ant} e^{ i \mathbf{k}(\nu, \phi, \theta) \cdot \mathbf{r}_{\rm ant}}
\mathcal{F}(\nu, \phi, \theta) = \sum_{\rm ant} w_{\rm ant} e^{ i \mathbf{k}(\nu, \phi, \theta) \cdot \mathbf{r}_{\rm ant}}
where :math:`\mathbf{k} = \frac{2\pi}{\lambda} (\cos \phi \cos \theta, \sin \phi \cos \theta, \sin \theta )`
is the wave vector for a wave propagation in a direction described by spherical coordinates,
:math:`\lambda` is the wavelength, :math:`\phi` is the azimuth,
:math:`\theta` is the elevation, :math:`\mathbf{r}_{\rm ant}`
is the antenna position matrix.
is the antenna position matrix, :math:`w_{\rm ant}` is the weight of the antenna (defined in :attr:`~nenupy.instru.interferometer.Interferometer.feed_weights`).
This method considers the ``sky`` as the desired output (in terms of
time, frequency and sky positions). It evaluates the effective
Expand All @@ -546,6 +590,12 @@ def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = Fals
:class:`~nenupy.astro.pointing.Pointing`
:param return_complex:
Return complex array factor if `True` or power if `False`
:type return_complex:
`bool`
:param normalize:
Return the normalized array factor. Default is `True`.
:type normalize:
`bool`
:return:
Array factor of the antenna distribution shaped as ``(time, frequency, 1, coordinates)``.
Expand All @@ -564,12 +614,7 @@ def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = Fals
# antenna position and the difference between sky and
# pointing ground projections.
geometric_delays = self._geometric_delays(sky, effective_pointing)

# Add delay errors for each antennae
# They could be cable connection errors during construction, cables of wrong length, ...
if self.extra_delays is not None:
geometric_delays += self.extra_delays[:, None, None, None, None]


# Use the sky frequency attribute to compute the wavelength
# and prepare the coefficient of the exponential with
# the correct dimensions.
Expand All @@ -579,30 +624,28 @@ def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = Fals
(1, 1, wavelength.size, 1, 1)
) # (antenna, time, frequency, polar, coord)

exponent = coeff * geometric_delays
exponent = coeff * (geometric_delays + self.antenna_delays[(...,) + (None,)*4])

# coord_chunk = exponent.shape[-1]//cpu_count()
# coord_chunk = 1 if coord_chunk == 0 else coord_chunk
# exponent = da.rechunk(
# exponent,
# chunks=exponent.shape[:-1] + (coord_chunk,)
# )
pre_sum = np.exp(exponent)
# apply fedd gain errors if any
if self.feed_gains is not None:
print(self.feed_gains)
pre_sum *= self.feed_gains[:, None, None, None, None]


complex_array_factor = np.sum(
pre_sum,
self.antenna_weights[(...,) + (None,)*4] * np.exp(exponent),
axis=0
)/exponent.shape[0] # normalized
)
if normalize:
complex_array_factor /= self.size

if return_complex:
return complex_array_factor
return complex_array_factor.real**2 + complex_array_factor.imag**2


def beam(self, sky: Sky, pointing: Pointing, return_complex: bool = False) -> Sky:
def beam(self, sky: Sky, pointing: Pointing, return_complex: bool = False, normalize: bool = True) -> Sky:
r""" Computes the phased-array response :math:`\mathcal{G}` over the ``sky`` for a given
``pointing``.
Expand Down Expand Up @@ -648,7 +691,8 @@ def beam(self, sky: Sky, pointing: Pointing, return_complex: bool = False) -> Sk
array_factor = self.array_factor(
sky=sky,
pointing=pointing,
return_complex=return_complex
return_complex=return_complex,
normalize=normalize
)

# Compute the total antenna gain, i.e. the sum of all
Expand Down
34 changes: 20 additions & 14 deletions nenupy/instru/nenufar.py
Original file line number Diff line number Diff line change
Expand Up @@ -288,7 +288,8 @@ class MiniArray(Interferometer):
~nenupy.instru.interferometer.Interferometer.antenna_gains
~nenupy.instru.interferometer.Interferometer.baselines
~nenupy.instru.interferometer.Interferometer.size
~nenupy.instru.interferometer.Interferometer.antenna_weights
~nenupy.instru.interferometer.Interferometer.antenna_delays
.. rubric:: Methods Summary
Expand All @@ -313,9 +314,10 @@ class MiniArray(Interferometer):
"""

def __init__(self,
index: int = 0,
extra_delays: np.ndarray = None,
feed_gains: np.ndarray = None, ):
index: int = 0,
antenna_delays: np.ndarray = None,
antenna_weights: np.ndarray = None
):
self.index = index

try:
Expand Down Expand Up @@ -350,12 +352,12 @@ def __init__(self,
])

super().__init__(position=position,
antenna_names=antenna_names,
antenna_positions=antenna_positions,
antenna_gains=antenna_gains,
extra_delays=extra_delays,
feed_gains=feed_gains
)
antenna_names=antenna_names,
antenna_positions=antenna_positions,
antenna_gains=antenna_gains,
antenna_delays=antenna_delays,
antenna_weights=antenna_weights
)


def __repr__(self):
Expand Down Expand Up @@ -412,7 +414,8 @@ def beam(self,
sky: Sky,
pointing: Pointing,
configuration: NenuFAR_Configuration = NenuFAR_Configuration(),
return_complex: bool = False
return_complex: bool = False,
normalize: bool = True
) -> Sky:
r""" Computes the Mini-Array beam over the ``sky`` for a given
``pointing``.
Expand Down Expand Up @@ -533,7 +536,8 @@ def beam(self,
return super().beam(
sky=sky,
pointing=self.analog_pointing(pointing, configuration=configuration),
return_complex=return_complex
return_complex=return_complex,
normalize=normalize
)# / aeff[None, :, None, None]


Expand Down Expand Up @@ -1142,7 +1146,8 @@ def beam(self,
pointing: Pointing,
analog_pointing: Pointing = None,
configuration: NenuFAR_Configuration = NenuFAR_Configuration(),
return_complex: bool = False
return_complex: bool = False,
normalize: bool = True
) -> Sky:
r""" Computes the NenuFAR beam over the ``sky`` for a given
``pointing``.
Expand Down Expand Up @@ -1239,7 +1244,8 @@ def beam(self,
sky=sky,
pointing=analog_pointing,
configuration=configuration,
return_complex=return_complex
return_complex=return_complex,
normalize=normalize
).value*count for gain, count in zip(self.antenna_gains[indices], counts)
]),
axis=0
Expand Down

0 comments on commit ef928ac

Please sign in to comment.