From fcc3d59ef28f7d39f5c0afce1e0c594f77aa5f69 Mon Sep 17 00:00:00 2001
From: Cedric Viou <cedric.viou@obs-nancay.fr>
Date: Wed, 22 Nov 2023 18:42:57 +0100
Subject: [PATCH] Add support for array imperfections

---
 nenupy/instru/interferometer.py | 33 ++++++++++++++++++++++++---------
 nenupy/instru/nenufar.py        | 18 +++++++++++-------
 2 files changed, 35 insertions(+), 16 deletions(-)

diff --git a/nenupy/instru/interferometer.py b/nenupy/instru/interferometer.py
index d53eaba..13f9c6d 100644
--- a/nenupy/instru/interferometer.py
+++ b/nenupy/instru/interferometer.py
@@ -205,15 +205,19 @@ class Interferometer(ABC, metaclass=CombinedMeta):
     """
 
     def __init__(self,
-            position: EarthLocation,
-            antenna_names: np.ndarray,
-            antenna_positions: np.ndarray,
-            antenna_gains: np.ndarray
-        ):
+                 position: EarthLocation,
+                 antenna_names: np.ndarray,
+                 antenna_positions: np.ndarray,
+                 antenna_gains: np.ndarray,
+                 extra_delays: np.ndarray = None,
+                 feed_gains: 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
 
 
     def __getitem__(self, n):
@@ -540,6 +544,8 @@ def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = Fals
                 the targeted pointing directions over the time.
             :type pointing:
                 :class:`~nenupy.astro.pointing.Pointing`
+            :param return_complex:
+                Return complex array factor if `True` or power if `False`
 
             :return:
                 Array factor of the antenna distribution shaped as ``(time, frequency, 1, coordinates)``.
@@ -558,7 +564,12 @@ 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.
@@ -569,6 +580,12 @@ def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = Fals
         ) # (antenna, time, frequency, polar, coord)
 
         exponent = coeff * geometric_delays
+        
+        # apply fedd gain errors if any
+        if self.feed_gains is not None:
+            print(self.feed_gains)
+            exponent *= self.feed_gains[:, None, None, None, None]
+            
         # coord_chunk = exponent.shape[-1]//cpu_count()
         # coord_chunk = 1 if coord_chunk == 0 else coord_chunk
         # exponent = da.rechunk(
@@ -582,9 +599,7 @@ def array_factor(self, sky: Sky, pointing: Pointing, return_complex: bool = Fals
 
         if return_complex:
             return complex_array_factor
-        else:
-            # return np.sqrt(np.real(complex_array_factor * np.conjugate(complex_array_factor)))
-            return np.real(complex_array_factor * np.conjugate(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:
diff --git a/nenupy/instru/nenufar.py b/nenupy/instru/nenufar.py
index cb33f1b..c97d161 100644
--- a/nenupy/instru/nenufar.py
+++ b/nenupy/instru/nenufar.py
@@ -312,7 +312,10 @@ class MiniArray(Interferometer):
 
     """
 
-    def __init__(self, index: int = 0):
+    def __init__(self,
+                 index: int = 0,
+                 extra_delays: np.ndarray = None,
+                 feed_gains: np.ndarray = None, ):
         self.index = index
 
         try:
@@ -346,12 +349,13 @@ def __init__(self, index: int = 0):
             self._antenna_gain for _ in range(antenna_names.size)
         ])
 
-        super().__init__(
-            position=position,
-            antenna_names=antenna_names,
-            antenna_positions=antenna_positions,
-            antenna_gains=antenna_gains
-        )
+        super().__init__(position=position,
+                         antenna_names=antenna_names,
+                         antenna_positions=antenna_positions,
+                         antenna_gains=antenna_gains,
+                         extra_delays=extra_delays,
+                         feed_gains=feed_gains
+                         )
 
 
     def __repr__(self):