Interpolated spectrum

data/config.ini

@@ -8,6 +8,7 @@ psize=Psize.dat
 raw_ts=cerema/raw/201702281700.raw,cerema/raw/201702281730.raw
 raw_spec=cerema/spt/201702281715.spt
 hires_step=0.5
+cycle=14400
 
 [out]
 margin=0.005

data/processing/spec.py

@@ -5,6 +5,7 @@ import pathlib
 
 import matplotlib.pyplot as plt
 import numpy as np
+from scipy.interpolate import griddata
 
 parser = argparse.ArgumentParser(description="Pre-process time-series")
 parser.add_argument("-v", "--verbose", action="count", default=0)
@@ -39,16 +40,24 @@ inp = np.loadtxt(
     max_rows=64,
 )
 
+cycle = config.getfloat("inp", "cycle", fallback=None)
+if cycle is None:
+    f = inp["f"]
+    S = inp["S"] * Sm
+else:
+    f = np.arange(inp["f"].min(), inp["f"].max() + 1/cycle, 1/cycle)
+    S = griddata(inp["f"], inp["S"] * Sm, f)
+
 with out_spec.open("w") as out:
     out.write("SPEC1D\n")
-    np.savetxt(out, np.stack((inp["f"], inp["S"] * Sm), axis=1))
+    np.savetxt(out, np.stack((f, S), axis=1))
 
-df = np.diff(inp['f']).min()
+df = np.diff(f).min()
 log.info(f"Minimum frequency delta: {df:.4f}Hz")
 log.info(f"Maximum modelled time: {1/df:.0f}s")
 
 fig, ax = plt.subplots()
-ax.plot(inp["f"], inp["S"]*Sm, c="k", lw=1)
+ax.plot(f, S, c="k", lw=1)
 ax.autoscale(True, "x", tight=True)
 ax.grid()
 ax.set(xlim=0, ylim=0, xlabel="f (Hz)", ylabel="S (m^2/Hz)")