LALPulsar  6.1.0.1-c9a8ef6
testSFTWindows.py
Go to the documentation of this file.
1 import sys
2 import csv
3 import subprocess
4 import numpy as np
5 
6 maxd = 2e-3
7 maxdrms = 1e-5
8 
9 exitcode = 0
10 
11 try:
12  import matplotlib as mpl
13 
14  mpl.use("Agg")
15  import matplotlib.pyplot as plt
16  import matplotlib.gridspec as gridspec
17 
18  plot = True
19 except:
20  plot = False
21 
22 # generate window data
23 window_data = subprocess.check_output(["./genSFTWindows"]).decode("utf-8").splitlines()
24 
25 # parse window data from CSV
26 window_csv = csv.reader(window_data)
27 window_csv_itr = iter(window_csv)
28 header = next(window_csv)
29 csv_cols = dict()
30 for name in header:
31  csv_cols[name] = list()
32 for row in window_csv_itr:
33  for name, value in zip(header, row):
34  csv_cols[name].append(float(value))
35 winrms = dict()
36 win = dict()
37 for name in csv_cols:
38  winrms[name] = csv_cols[name][0]
39  win[name] = np.array(csv_cols[name][1:])
40 t = np.arange(len(win[name])) / 256
41 
42 # compare windows
43 for k in ("Tukey", "Hann"):
44  win_ref_name, win_ref = [
45  (name, win[name]) for name in win if ("Create" + k) in name
46  ][0]
47  if plot:
48  headtail = 512
49  ranges = [
50  range(0, len(t)),
51  range(0, headtail),
52  range(len(t) - headtail, len(t)),
53  ]
54  legend_loc = ["lower center", "lower right", "lower left"]
55  fig = plt.figure(tight_layout=True, figsize=(16, 8))
56  gs = gridspec.GridSpec(3, 2)
57  axs = [
58  fig.add_subplot(gs[0, :]),
59  fig.add_subplot(gs[1, 0]),
60  fig.add_subplot(gs[1, 1]),
61  fig.add_subplot(gs[2, 0]),
62  fig.add_subplot(gs[2, 1]),
63  ]
64  for i, ax in enumerate(axs[0:3]):
65  ax.plot(
66  t[ranges[i]],
67  win_ref[ranges[i]],
68  label=win_ref_name,
69  color="black",
70  linewidth=2,
71  zorder=-10,
72  )
73  ax.set_xticks(np.linspace(min(t[ranges[i]]), max(t[ranges[i]] + t[1]), 5))
74  ax.ticklabel_format(axis="x", useOffset=False)
75  ax.ticklabel_format(axis="y", style="sci", scilimits=(0, 0))
76  ax.set_xlabel("Time / seconds")
77  ax.set_ylabel("Window value")
78  for n, (win_cmp_name, win_cmp) in enumerate(
79  [(name, win[name]) for name in win if (k in name and name != win_ref_name)]
80  ):
81  if plot:
82  for i, ax in enumerate(axs[0:3]):
83  ax.plot(
84  t[ranges[i]], win_cmp[ranges[i]], label=win_cmp_name, linewidth=1
85  )
86  ax.legend(loc=legend_loc[i])
87  ax.ticklabel_format(axis="x", useOffset=False)
88  ax.ticklabel_format(axis="y", style="sci", scilimits=(0, 0))
89  d = win_cmp - win_ref
90  md = max(abs(d))
91  print(f"max(|{win_cmp_name} - {win_ref_name}|) = {md:0.2e}, max = {maxd:0.2e}")
92  if md > maxd:
93  exitcode = 1
94  drms = abs(winrms[win_cmp_name] - winrms[win_ref_name])
95  print(
96  f"abs(|{win_cmp_name}.rms - {win_ref_name}.rms|) = {drms:0.2e}, max = {maxdrms:0.2e}"
97  )
98  if drms > maxdrms:
99  exitcode = 1
100  if plot:
101  for i, ax in enumerate(axs[3:]):
102  ax.plot(t[ranges[i + 1]], d[ranges[i + 1]], color="black")
103  ax.ticklabel_format(axis="x", useOffset=False)
104  ax.ticklabel_format(axis="y", style="sci", scilimits=(0, 0))
105  ax.set_xticks(
106  np.linspace(min(t[ranges[i + 1]]), max(t[ranges[i + 1]] + t[1]), 5)
107  )
108  ax.ticklabel_format(axis="x", useOffset=False)
109  ax.ticklabel_format(axis="y", style="sci", scilimits=(0, 0))
110  ax.set_xlabel("Time / seconds")
111  ax.set_ylabel("Window difference")
112  plt.savefig(f"SFTWindows{k}{n}.png")
113 
114 sys.exit(exitcode)
#define min(a, b)
#define max(a, b)