|
| 1 | +################################################################################ |
| 2 | +# Copyright (c) 2026, National Research Foundation (SARAO) |
| 3 | +# |
| 4 | +# Licensed under the BSD 3-Clause License (the "License"); you may not use |
| 5 | +# this file except in compliance with the License. You may obtain a copy |
| 6 | +# of the License at |
| 7 | +# |
| 8 | +# https://opensource.org/licenses/BSD-3-Clause |
| 9 | +# |
| 10 | +# Unless required by applicable law or agreed to in writing, software |
| 11 | +# distributed under the License is distributed on an "AS IS" BASIS, |
| 12 | +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. |
| 13 | +# See the License for the specific language governing permissions and |
| 14 | +# limitations under the License. |
| 15 | +################################################################################ |
| 16 | + |
| 17 | +"""Sample test for tied-array-resampled-voltage stream.""" |
| 18 | + |
| 19 | +import asyncio |
| 20 | +from collections.abc import AsyncGenerator |
| 21 | + |
| 22 | +import aiokatcp |
| 23 | +import numpy as np |
| 24 | +import pytest |
| 25 | +from pytest_check import check |
| 26 | + |
| 27 | +from katgpucbf.pytest_plugins.reporter import Reporter |
| 28 | +from katgpucbf.utils import TimeConverter |
| 29 | +from qualification.cbf import CBFRemoteControl |
| 30 | + |
| 31 | +from ..recv import TiedArrayChannelisedVoltageReceiver, TiedArrayResampledVoltageReceiver |
| 32 | + |
| 33 | + |
| 34 | +@pytest.fixture |
| 35 | +async def sensor_watcher(cbf: CBFRemoteControl) -> AsyncGenerator[aiokatcp.SensorWatcher, None]: |
| 36 | + """Establish a secondary connection to the product controller with a sensor watcher. |
| 37 | +
|
| 38 | + The yielded sensor watcher is not yet synchronised. |
| 39 | + """ |
| 40 | + # aiokatcp doesn't currently handle adding watchers after the connection |
| 41 | + # is already established; SensorWatcher is also somewhat expensive. So |
| 42 | + # instead we create a separate connection for monitoring sensors. |
| 43 | + secondary = aiokatcp.Client(*cbf.product_controller_endpoint) |
| 44 | + sensor_watcher = aiokatcp.SensorWatcher(secondary) |
| 45 | + secondary.add_sensor_watcher(sensor_watcher) |
| 46 | + |
| 47 | + yield sensor_watcher |
| 48 | + |
| 49 | + secondary.close() |
| 50 | + await secondary.wait_closed() |
| 51 | + |
| 52 | + |
| 53 | +@pytest.mark.name("VLBI mean power") |
| 54 | +async def test_mean_power( |
| 55 | + pdf_report: Reporter, |
| 56 | + receive_tied_array_resampled_voltage: TiedArrayResampledVoltageReceiver | None, |
| 57 | + receive_tied_array_channelised_voltage: TiedArrayChannelisedVoltageReceiver, |
| 58 | + cbf: CBFRemoteControl, |
| 59 | + sensor_watcher: aiokatcp.SensorWatcher, |
| 60 | + pass_channels: slice, |
| 61 | +) -> None: |
| 62 | + """Test mean-power sensor values against tied-array channelised voltage. |
| 63 | +
|
| 64 | + Verification method |
| 65 | + ------------------- |
| 66 | + Verified by means of test. Inject a white noise signal and set beam weights |
| 67 | + so that a single antenna contributes. Wait until each ``mean-power`` sensor |
| 68 | + timestamp is after the system steady-state timestamp (plus one |
| 69 | + ``power-int-time`` so the averaging window is entirely post-steady-state). |
| 70 | + Measure mean power from the tied-array channelised voltage stream over the |
| 71 | + passband channels, and compare against each ``mean-power`` sensor. |
| 72 | + """ |
| 73 | + assert receive_tied_array_resampled_voltage is not None |
| 74 | + receiver = receive_tied_array_resampled_voltage |
| 75 | + pcc = cbf.product_controller_client |
| 76 | + |
| 77 | + pdf_report.step("Setup signal generator and gains.") |
| 78 | + async with asyncio.TaskGroup() as tg: |
| 79 | + for i, name in enumerate(receive_tied_array_channelised_voltage.stream_names): |
| 80 | + gains = [0.0] * len(receive_tied_array_channelised_voltage.source_indices[i]) |
| 81 | + gains[0] = 1.0 |
| 82 | + tg.create_task(pcc.request("beam-weights", name, *gains)) |
| 83 | + pdf_report.detail(f"Set beam weights to {gains}.") |
| 84 | + |
| 85 | + async with asyncio.TaskGroup() as tg: |
| 86 | + for dsim_name in cbf.dsim_names: |
| 87 | + tg.create_task(pcc.request("dsim-signals", dsim_name, "common=wgn(0.02);common;common;")) |
| 88 | + pdf_report.detail("Set dsim signals white noise.") |
| 89 | + |
| 90 | + pdf_report.step("Wait for mean-power sensors to reach steady state.") |
| 91 | + await sensor_watcher.synced.wait() # Implicitly waits for connection too |
| 92 | + time_converter = TimeConverter(receiver.sync_time, receiver.scale_factor_timestamp) |
| 93 | + steady_state_unix = time_converter.adc_to_unix(await cbf.steady_state_timestamp()) |
| 94 | + # Require a full power-int-time of data after steady state so the sensor |
| 95 | + # average does not include pre-change samples. |
| 96 | + min_sensor_time = steady_state_unix + receiver.power_int_time |
| 97 | + |
| 98 | + sensor_names = [ |
| 99 | + f"{receiver.stream_names[0]}.{pol}{chan}.mean-power" |
| 100 | + for pol in receiver.pol_ordering |
| 101 | + for chan in range(receiver.n_chans) |
| 102 | + ] |
| 103 | + |
| 104 | + async def wait_mean_power_steady_state() -> None: |
| 105 | + while True: |
| 106 | + timestamps = [sensor_watcher.sensors[name].timestamp for name in sensor_names] |
| 107 | + earliest = min(timestamps) |
| 108 | + if earliest >= min_sensor_time: |
| 109 | + pdf_report.detail("Mean-power sensors reached steady state timestamp.") |
| 110 | + break |
| 111 | + await asyncio.sleep(0.5) |
| 112 | + |
| 113 | + await asyncio.wait_for(asyncio.create_task(wait_mean_power_steady_state()), timeout=15.0) |
| 114 | + |
| 115 | + pdf_report.step("Measure power from tied-array channelised voltage.") |
| 116 | + _, tacv_data = await receive_tied_array_channelised_voltage.next_complete_chunk() |
| 117 | + tacv_data = tacv_data.astype(np.float64).view(np.complex128)[..., 0] # Convert to complex128 |
| 118 | + # Only use the pass channels for beam zero for the power calculation. |
| 119 | + tacv_data = tacv_data[0][pass_channels] |
| 120 | + tacv_power = (np.square(tacv_data.real) + np.square(tacv_data.imag)).mean() |
| 121 | + pdf_report.detail(f"Mean TACV power over passband channels: {tacv_power}.") |
| 122 | + # Test that we aren't accidentally testing zero values: |
| 123 | + assert tacv_power > 0.0 |
| 124 | + |
| 125 | + pdf_report.step("Compare mean-power sensors against TACV power.") |
| 126 | + for sensor_name in sensor_names: |
| 127 | + sensor = sensor_watcher.sensors[sensor_name] |
| 128 | + with check: |
| 129 | + assert sensor.timestamp >= min_sensor_time |
| 130 | + assert sensor.value == pytest.approx(tacv_power, rel=5e-3), ( |
| 131 | + f"TACV power ^2: {tacv_power} does not match total theta^2: {sensor.value}" |
| 132 | + + f" for sensor {sensor_name}" |
| 133 | + ) |
| 134 | + pdf_report.detail("Power agrees to within 0.5%.") |
0 commit comments