From 52c61ee9f6cbb38d119f7878d9f96e191ea9d168 Mon Sep 17 00:00:00 2001
From: Eric Kooistra <kooistra@astron.nl>
Date: Wed, 11 Dec 2024 13:11:03 +0100
Subject: [PATCH] Update TODO in non_maximal_upsample_bpf().
---
.../lofar2/model/rtdsp/singlechannel.py | 102 ++++++++----------
1 file changed, 46 insertions(+), 56 deletions(-)
diff --git a/applications/lofar2/model/rtdsp/singlechannel.py b/applications/lofar2/model/rtdsp/singlechannel.py
index ced9bcf3b2..e51cbac059 100644
--- a/applications/lofar2/model/rtdsp/singlechannel.py
+++ b/applications/lofar2/model/rtdsp/singlechannel.py
@@ -264,28 +264,52 @@ def non_maximal_downsample_bpf(x, Ndown, k, Ndft, coefs, verbosity=1):
def non_maximal_upsample_bpf(xBase, Nup, k, Ndft, coefs, verbosity=1):
"""BPF xBase at bin k in range(Nup), and upsample xBase by factor U = Nup.
- Implement nonmaximal upsampling upconverter for one bin. The maximal
- upsampling upconverter for one bin has the kPhasors per polyphase
- branch similar of [HARRIS Fig 7.16]. For nonmaximal this needs to be
- extended with the tPhasors per polyphase branch similar as for down in
- [HARRIS Eq. 9.3], to compensate for the oversampling time shift.
-
- The BPF is centered at w_k = 2pi k / Ndft, where Ndft is number of
- frequency bins, is DFT size. The polyphase FIR structure has Nphases = Ndft
- branches, to fit the requested number of bins. The polyphase FIR structure
- is maximally upsampled (= critically sampled) for Nup = Ndft, but it
- can support any Nup <= Ndft. The output data shifts out per Nup samples,
- so it appears from different branches when Nup < Ndft and a new block is
- shifted out. Therefore the output data cannot be FIR filtered per branch
- for the whole input xBase. Instead it needs to be FIR filtered per block of
- Nup output samples, using pfs.polyDelays in pfs.filter_up().
-
TODO
- . This code only runs ok for Ros = 1, 2 when Ndft = 16, but not for
- Ros = 4, so need to fix that. According to [CROCHIERE 7.2.7] the
- polyphase structure is only suitable for Ros is integer >= 1. For other
- Ros > 1 the weighted overlap-add (WOLA) structure is suitable, so
- need to add WOLA.
+ . According to [CROCHIERE 7.2.7] the polyphase structure is only suitable
+ for Ros is integer >= 1.
+ . For other Ros > 1 the weighted overlap-add (WOLA) structure is suitable,
+ so need to add WOLA.
+
+ Remark:
+ . It is difficult to transpose the WOLA approach to a PFS approach that
+ suits any Ros > 1. For example, with blocks a, b, c, d, f, ... of Ndft
+ samples and block a = a0 a1 for Ros = 2, and replicating the block Ntaps
+ = 2 times in time, WOLA using rows:
+
+ a0 a1 a0 a1
+ b0 b1 b0 b1
+ c0 c1 c0 c1
+ d0 d1 d0 d1
+ e0 e1 e0 e1
+ f0 f1 f0 f1
+ g0 g1 g0 g1 --> time
+
+ results in PFS using columns:
+
+ d0 e0 f0 g0
+ c1 d1 e1 f1
+ b0 c0 d0 e0
+ a1 b1 c1 d1
+
+ For fractional Ros, e.g. Ros = 4/3, it can be tried to find a pattern
+ processing per column, but that is much more complicated then processing
+ per row as with WOLA (Ntaps = 2):
+
+ a0 a1 a2 a3 a0 a1 a2 a3
+ b0 b1 b2 b3 b0 b1 b2 b3
+ c0 c1 c2 c3 c0 c1 c2 c3
+ d0 d1 d2 d3 d0 d1 d2 d3
+ e0 e1 e2 e3 e0 e1 e2 e3
+ f0 f1 f2 f3 f0 f1 f2 f3
+ results in PFS using columns:
+
+ c0 c1 c2 d0 d1 d2 e0 e1 e2 f0 f1
+ b3 b0 b1 c3 c0 c1 d3 d0 d1 e3 e0
+ a2 a3 b2 b3 c2 c3 d2 d3
+
+ Extra complication is for e.g. Ros = 32 / 27, because 32 does not divide
+ (32 - 27). Anyway, for fractional Ros PFS approach seems not useful,
+ because the PFS is then changing for every block in time.
Input:
. xBase: Input equivalent baseband signal xBase[m]
@@ -305,46 +329,12 @@ def non_maximal_upsample_bpf(xBase, Nup, k, Ndft, coefs, verbosity=1):
# print('WARNING: Only support integer Ros')
# return None
- Nblocks = len(xBase)
-
- # Prepare output
- polyYc = np.zeros((Nup, Nblocks), dtype='cfloat')
-
- # Adjust LO phase for group delay of analysis LPF and synthesis LPF in
- # series. This is not needed in practise, but is needed to be able to
- # exactly compare reconstructed signal with original input time series
- # signal. Assume coefs has group delay (len(coefs) - 1) / 2.
- hPairGroupDelay = len(coefs) - 1
- hPhasor = np.exp(-2j * np.pi * k / Ndft * hPairGroupDelay)
-
- # Oversampling time shift compensation via frequency dependent phase shift
- tPhasors = time_shift_phasor_arr(k, Nup, Ndft, Nblocks, -1) # [HARRIS Eq 9.3]
-
- # PFS with Ndft polyphases
- pfs = PolyPhaseFirFilterStructure(Ndft, coefs, cmplx=True)
- kPhasors = unit_circle_loops_phasor_arr(k, Ndft, 1) # [HARRIS Eq 7.8]
-
- # Polyphase FIR filter input xBase
- for b in range(Nblocks):
- # Phase rotate polyphases for bin k
- inData = xBase[b] * hPhasor * tPhasors[b] * kPhasors
-
- # Filter block [HARRIS Fig 9.15, CROCHIERE Fig 7.16]]
- pfsData = pfs.filter_block(inData)
-
- # Output Nup samples yc[n] for every input sample xBase[m]
- polyYc[:, b] = Ndft * pfsData
-
- yc = polyYc.T.reshape(1, Nup * Nblocks)[0]
-
if verbosity:
print('> non_maximal_upsample_bpf():')
print(' . len(xBase) =', str(len(xBase)))
- print(' . Nblocks =', str(Nblocks))
- print(' . len(yc) =', str(len(yc))) # = Nblocks
print(' . Ros =', str(Ros))
print(' . Nup =', str(Nup))
print(' . Ndft =', str(Ndft))
print(' . k =', str(k))
print('')
- return yc
+ return None
--
GitLab