irfft: apply ramp function to prevent double-counting in asymmetric ifgs

Author: stuart-cls

orangecontrib/spectroscopy/irfft.py

@@ -55,6 +55,12 @@ def find_zpd(ifg, peak_search):
     else:
         raise NotImplementedError
 
+def wing_size(ifg, zpd):
+    """Calculate negative and positive wing size (including zpd point)"""
+    wing_n = zpd + 1
+    wing_p = ifg.shape[-1] - zpd
+    return wing_n, wing_p
+
 def apodize(ifg, zpd, apod_func):
     """
     Perform apodization of asymmetric interferogram using selected apodization
@@ -80,8 +86,7 @@ def apodize(ifg, zpd, apod_func):
     # Calculate negative and positive wing size
     # correcting zpd from 0-based index
     ifg_N = ifg.shape[-1]
-    wing_n = zpd + 1
-    wing_p = ifg_N - zpd
+    wing_n, wing_p = wing_size(ifg, zpd)
     # Use larger wing to define apodization window width
     M = max(wing_n, wing_p)
 
@@ -132,6 +137,18 @@ def apodize(ifg, zpd, apod_func):
 
     return ifg_apod
 
+def ramp(ifg, zpd):
+    wing_n, wing_p = wing_size(ifg, zpd)
+    # Use smaller wing to define symmetric subset size
+    M = min(wing_n, wing_p)
+    ramp = np.linspace(0, 1, 2 * M - 1)
+    if wing_n < wing_p:
+        ifg[..., :ramp.shape[0]] *= ramp
+    else:
+        ifg[..., -ramp.shape[0]:] *= ramp[::-1]
+
+    return ifg
+
 def _zero_fill_size(ifg_N, zff):
     # Calculate desired array size
     Nzff = ifg_N * zff
@@ -215,6 +232,7 @@ def __call__(self, ifg, zpd=None, phase=None):
             # Note that L is now the zpd index
             Ixs = ifg[self.zpd - L : self.zpd + L].copy()
             ifg = apodize(ifg, self.zpd, self.apod_func)
+            ifg = ramp(ifg, self.zpd)
             ifg = zero_fill(ifg, self.zff)
             ifg = np.hstack((ifg[self.zpd:], ifg[0:self.zpd]))
             Nzff = ifg.shape[0]
@@ -300,6 +318,7 @@ def __call__(self, ifg, zpd=None, phase=None):
             # Note that L is now the zpd index
             Ixs = ifg[:, self.zpd - L : self.zpd + L].copy()
             ifg = apodize(ifg, self.zpd, self.apod_func)
+            ifg = ramp(ifg, self.zpd)
             ifg = zero_fill(ifg, self.zff)
             ifg = np.hstack((ifg[:, self.zpd:], ifg[:, 0:self.zpd]))
             Nzff = ifg.shape[1]

orangecontrib/spectroscopy/tests/test_irfft.py

@@ -7,7 +7,7 @@
 
 from orangecontrib.spectroscopy.irfft import (IRFFT, zero_fill, PhaseCorrection,
                                               find_zpd, PeakSearch, ApodFunc,
-                                              MultiIRFFT, apodize,
+                                              MultiIRFFT, apodize, ramp,
                                               )
 
 dx = 1.0 / 15797.337544 / 2.0
@@ -19,6 +19,14 @@ def setUp(self):
         self.ifg_single = Orange.data.Table("IFG_single.dpt")
         self.ifg_seq_ref = Orange.data.Table("agilent/background_agg256.seq")
         self.sc_dat_ref = Orange.data.Table("agilent/background_agg256.dat")
+        self.ifgs = {
+            'even_sym': self.ifg_single.X[0],
+            'odd_sym': self.ifg_single.X[0][1:],
+            'even_asym': self.ifg_seq_ref.X[0][1:],
+            'odd_asym': self.ifg_seq_ref.X[0],
+            'even_asym_reverse': self.ifg_seq_ref.X[0][1:][::-1],
+            'odd_asym_reverse': self.ifg_seq_ref.X[0][::-1],
+        }
 
     def test_zero_fill(self):
         N = 1975
@@ -144,18 +152,19 @@ def test_multi_ab(self):
         np.testing.assert_allclose(ab[:, limits[0]:limits[1]], dat, atol=0.004)
 
     def test_apodization(self):
-        ifgs = {
-            'even_sym': self.ifg_single.X[0],
-            'odd_sym': self.ifg_single.X[0][1:],
-            'even_asym': self.ifg_seq_ref.X[0][1:],
-            'odd_asym': self.ifg_seq_ref.X[0],
-            'even_asym_reverse': self.ifg_seq_ref.X[0][1:][::-1],
-            'odd_asym_reverse': self.ifg_seq_ref.X[0][::-1],
-        }
         for apod_func in ApodFunc:
-            for k, ifg in ifgs.items():
+            for k, ifg in self.ifgs.items():
                 with self.subTest(apod_func=apod_func.name, ifg=k):
                     zpd = find_zpd(ifg, PeakSearch.ABSOLUTE)
                     out = apodize(ifg, zpd, apod_func)
                     # Apodization should not change value at zpd
                     self.assertAlmostEqual(ifg[zpd], out[zpd])
+
+    def test_ramp(self):
+        ifg = self.ifgs['odd_asym']
+        zpd = find_zpd(ifg, PeakSearch.ABSOLUTE)
+        ramp_fwd = ramp(ifg, zpd)
+        ifg_rev = self.ifgs['odd_asym_reverse']
+        zpd_rev = find_zpd(ifg_rev, PeakSearch.ABSOLUTE)
+        ramp_rev = ramp(ifg_rev, zpd_rev)
+        np.testing.assert_array_equal(ramp_fwd, ramp_rev[::-1])