Detect noise level

zonnewaarneming.py

@@ -34,8 +34,26 @@ def set_stoprequest(event):
     global stoprequest
     stoprequest = True
 
+def update_logticks(ax, noiselevel):
+    yticks = []
+    yticklabels = []
+    ytick = noiselevel
+    db = 0
+    while ytick < ax.get_ylim()[1]:
+        yticks.append(ytick)
+        yticklabels.append(f"{db} dB")
+        ytick *= 10**0.1
+        db += 1
+    while len(yticks) > 15:
+        yticks = yticks[::2]
+        yticklabels = yticklabels[::2]
+    ax.set_yticks(yticks)
+    ax.set_yticklabels(yticklabels)
+
 signal = 0.
 prev_signal = 0.
+noise_sum = 0.
+noise_count = 0
 
 y_scale_is_logarithmic = False
 
@@ -164,6 +182,10 @@ while not stoprequest:
         if separation_alt < 0*u.deg:
             separation = -separation
 
+    if dt.el > 20*u.deg and abs(separation) > 2*u.deg and closest_source != "Zon":
+        noise_sum += signal
+        noise_count += 1
+
     if abs(dt.dist_el) > 0.2 * u.deg:
         # Only show values if the elevation is already on source
         print(f"Elevation off ({separation_alt}), continuing")
@@ -182,7 +204,11 @@ while not stoprequest:
         maxsignal = signal
     if signal < minsignal:
         minsignal = signal
-    noiselevel = 75 # TODO: infer this from data
+    if noise_count < 10:
+        noiselevel = 75 # TODO: infer this from data
+    else:
+        noiselevel = noise_sum / noise_count
+        print(f"New noise level: {noiselevel}")
 
     scatter_data = scatterplot.get_offsets()
     scatter_data = np.vstack((scatter_data, [separation.deg, signal]))
@@ -191,29 +217,17 @@ while not stoprequest:
     prev_points = ax.scatter([separation.deg], [signal], s=50, color='red')
     ax.set_ylim((0.9 * minsignal, 1.1 * maxsignal))
 
-    if (not y_scale_is_logarithmic) and minsignal < 1000:
+    if minsignal < 1000:
         # Attenuator on
         y_scale_is_logarithmic = True
         ax.set_yscale("log")
 
-        yticks = []
-        yticklabels = []
-        ytick = noiselevel
-        db = 0
-        while ytick < ax.get_ylim()[1]:
-            yticks.append(ytick)
-            yticklabels.append(f"{db} dB")
-            ytick *= 10**0.1
-            db += 1
-        while len(yticks) > 15:
-            yticks = yticks[::2]
-            yticklabels = yticklabels[::2]
-        ax.set_yticks(yticks)
-        ax.set_yticklabels(yticklabels)
         ax.get_yaxis().set_tick_params(which='minor', size=0)
         ax.get_yaxis().set_tick_params(which='minor', width=0)
         ax.set_ylabel("Relatieve intensiteit")
 
+        update_logticks(ax, noiselevel)
+
     update_legend(ax, closest_source)
 
     # show the plot