#!/usr/bin/env python3
import sys
import numpy as np
if len(sys.argv) != 2:
print("Must supply the freq resolution in kHz as first argument (Must be a valid MWAX fine channel width)")
Fs = float(sys.argv[1]) #Fine Channel Resolution in kHz
if (Fs*1000) % 200 != 0:
print("Fine channel resolution must be a multiple of 200Hz")
exit()
if 1280000 % (Fs*1000) != 0:
print("Fine channel resolution must divide into 1280 kHz (Integer number of fine-channels per coarse channel)")
exit()
# --- Compute Frequency Response (at 200Hz ultra-fine channels) ---
h = np.loadtxt('MWARX_RRI_PrototypeFilter_512x8.dat')
c = int(327680000/200) # DC bin index (centre bin)
o = int(640000/200) # Size of half a channel
H = np.fft.fftshift(np.fft.fft(h,2*c))[c-4*o:c+4*o] # Do FFT & slice +-2.56MHz of data
# --- Compute Aliased Frequency Response ---
Hp = np.abs(H)**2 # Compute Power Response
Ha = Hp + np.roll(Hp,2*o) + np.roll(Hp,-2*o) # Add aliased components by shifting response by 1.28MHz in pos & neg directions.
Ha = Ha[3*o:-3*o] # Slice +-0.64MHz of data, i.e. A single 1.28MHz wide channel.
# --- MWAX Fscrunch ---
fscrunch = int(Fs/0.2)
fcentre = int(fscrunch/2)
fbins = int(1280/Fs)
tmp = np.reshape( np.roll(Ha,fcentre),(fbins,fscrunch)) # Group ultra-fine channels into fine-channel groups
if fscrunch % 2 == 0:
tmp[:,0] = ( tmp[:,0] + np.roll(Ha,-fcentre)[0::fscrunch] ) / 2 # If fscrunch is even: Combine edge channels together using an average
Ha = np.sum(tmp,axis=1) # Sum ultra-fine channels to get fine-channels
Ha = Ha / np.abs(Ha[int(fbins/2)]) # Normalise Magnitude (based on DC bin)
# --- Save to stdout as txt ---
np.savetxt( sys.stdout.buffer, Ha)