Merge pull request #1221 from qiboteam/set-if-tof

Adding detuning to `time_of_flight` experiment

Author: andrea-pasquale

doc/source/protocols/signal/time_of_flight.png

@@ -24,42 +24,27 @@ How to execute the experiment
     - id: time of flight experiment      # custom ID of the experiment
       operation: time_of_flight_readout  # unique name of the routine
       parameters:
-        readout_amplitude: 1             # usually high
+        readout_amplitude: 0.5
         window_size: 10
+        detuning: 50_000_000
         nshots: 1024
         relaxation_time: 20_000
 
 Although it is possible to avoid setting a specific readout amplitude, it is generally useful to set a high value here. Indeed, we are not looking for the optimal amplitude but we want to have a signal with enough power so that it is clear when it starts.
 
+
 Acquisition
 ^^^^^^^^^^^
 
-.. testcode::
-   :hide:
-
-   from qibolab import AcquisitionType
-
-It is important to note that this experiment makes use of the RAW acquisition mode (see `Qibolab documentation <https://qibo.science/qibolab/stable/api-reference/qibolab.html#qibolab.execution_parameters.AcquisitionType>`_), which may require some specific care depending on the instrument employed (for some devices demodulation could be used, or this mode could be available for just a single qubit at a time).
-
-For a 3D cavity we expect a plot (considering demodulation) as the following:
-
 .. image:: time_of_flight.png
 
-In this image, for example, we can see that the fit procedure did not find the correct pulse starting point (that is around 200). This could be improved increasing the window size.
-
-For a 2D resonator, some more oscillation could appear:
-
-.. image:: time_of_flight_2d.png
-
-Just as an example, we provide here also a plot for a 2D resonator case, where the demodulation is not carried out.
-As we can see the plot is very different, but the starting point of the signal is still clearly visible.
-
-.. image:: time_of_flight_2d_nodem.png
-
 Fit
 ^^^
 
-The fit procedure (:func:`qibocal.protocols.signal_experiments.time_of_flight_readout._fit`) employs a moving average, returning the time when it is maximum, namely when the signal starts being acquired.
+The fit procedure computes the expected time at which the signal should appear.
+To estimate the time of flight in the fitting a threshold is estimated to distinguish the noise from the signal, then
+the first point where the signal exceed this value is selected as the time of flight.
+
 
 Requirements
 ^^^^^^^^^^^^

doc/source/protocols/signal/time_of_flight.rst

@@ -11,6 +11,8 @@
 from qibocal.auto.operation import Data, Parameters, QubitId, Results, Routine
 from qibocal.calibration import CalibrationPlatform
 
+from .utils import _get_lo_frequency
+
 __all__ = ["calibrate_state_discrimination"]
 
 
@@ -34,16 +36,6 @@ class CalibrateStateDiscriminationParameters(Parameters):
 """Custom dtype for CalibrateStateDiscrimination."""
 
 
-def _get_lo_frequency(platform: CalibrationPlatform, qubit: QubitId) -> float:
-    """Get LO frequency given QubitId.
-
-    Currently it assumes that instruments with LOs is first one.
-    """
-    probe = platform.channels[platform.qubits[qubit].probe]
-    lo_config = platform.config(probe.lo)
-    return lo_config.frequency
-
-
 @dataclass
 class CalibrateStateDiscriminationResults(Results):
     """Calibrate State Discrimination outputs."""

doc/source/protocols/signal/time_of_flight_2d.png

@@ -12,6 +12,8 @@
 from qibocal.result import magnitude
 from qibocal.update import replace
 
+from .utils import _get_lo_frequency
+
 __all__ = ["time_of_flight_readout"]
 
 MINIMUM_TOF = 24
@@ -22,6 +24,8 @@
 class TimeOfFlightReadoutParameters(Parameters):
     """TimeOfFlightReadout runcard inputs."""
 
+    detuning: float = 10e6
+    """Detuning with respect to corresponding LO frequency [Hz]."""
     readout_amplitude: Optional[int] = None
     """Amplitude of the readout pulse."""
     window_size: Optional[int] = 10
@@ -36,16 +40,14 @@ class TimeOfFlightReadoutResults(Results):
     """Raw fitting output."""
 
 
-TimeOfFlightReadoutType = np.dtype([("samples", np.float64)])
-
-
 @dataclass
 class TimeOfFlightReadoutData(Data):
     """TimeOfFlightReadout acquisition outputs."""
 
     windows_size: int
     sampling_rate: int
-
+    intermediate_frequency: float
+    amplitude: dict[QubitId, float] = field(default_factory=dict)
     data: dict[QubitId, npt.NDArray] = field(default_factory=dict)
     """Raw data acquired."""
 
@@ -57,34 +59,47 @@ def _acquisition(
 ) -> TimeOfFlightReadoutData:
     """Data acquisition for time of flight experiment."""
     sequence = PulseSequence()
-    ro_pulses = {}
     native = platform.natives.single_qubit
-    ro_channels = []
     for qubit in targets:
         ro_channel, ro_pulse = native[qubit].MZ()[0]
-        ro_channels.append(ro_channel)
         if params.readout_amplitude is not None:
-            ro_pulse = replace(ro_pulse, amplitude=params.readout_amplitude)
-        ro_pulses[qubit] = ro_pulse
+            probe = replace(ro_pulse.probe, amplitude=params.readout_amplitude)
+            ro_pulse = replace(ro_pulse, probe=probe)
         sequence.append((ro_channel, ro_pulse))
+
     results = platform.execute(
         [sequence],
         nshots=params.nshots,
         relaxation_time=params.relaxation_time,
         acquisition_type=AcquisitionType.RAW,
         averaging_mode=AveragingMode.CYCLIC,
-        updates=[{ro_channel: {"delay": MINIMUM_TOF} for ro_channel in ro_channels}],
+        updates=[
+            {
+                platform.qubits[qubit].acquisition: {"delay": MINIMUM_TOF},
+                platform.qubits[qubit].probe: {
+                    "frequency": _get_lo_frequency(platform, qubit) + params.detuning,
+                },
+            }
+            for qubit in targets
+        ],
     )
 
     data = TimeOfFlightReadoutData(
-        windows_size=params.window_size, sampling_rate=platform.sampling_rate
+        windows_size=params.window_size,
+        sampling_rate=platform.sampling_rate,
+        intermediate_frequency=params.detuning,
+        amplitude={
+            qubit: list(sequence.channel(platform.qubits[qubit].acquisition))[
+                -1
+            ].probe.amplitude
+            for qubit in targets
+        },
     )
-
     # retrieve and store the results for every qubit
     for qubit in targets:
-        samples = magnitude(results[ro_pulses[qubit].id])
-        # store the results
-        data.register_qubit(TimeOfFlightReadoutType, (qubit), dict(samples=samples))
+        acq_handle = list(sequence.channel(platform.qubits[qubit].acquisition))[-1].id
+        data.data[qubit] = results[acq_handle]
+
     return data
 
 
@@ -96,16 +111,27 @@ def _fit(data: TimeOfFlightReadoutData) -> TimeOfFlightReadoutResults:
 
     window_size = data.windows_size
     sampling_rate = data.sampling_rate
-
     for qubit in qubits:
-        qubit_data = data[qubit]
-        samples = qubit_data.samples
-        window_size = int(len(qubit_data) / 10)
-        th = (np.mean(samples[:window_size]) + np.mean(samples[:-window_size])) / 2
-        delay = np.where(samples > th)[0][0]
-        time_of_flight_readout = float(delay / sampling_rate + MINIMUM_TOF)
+        delays = []
+        for i in range(2):
+            samples = data.data[qubit][:, i]
+            window_size = int(len(samples) / 10)
+
+            for feat in ["min", "max"]:
+                th = (
+                    getattr(np, feat)(samples[:window_size])
+                    + getattr(np, feat)(samples[:-window_size])
+                ) / 2
+                # try-expect in order to handle sporadic failing with mock
+                try:
+                    delay = np.where(samples < th if feat == "min" else samples > th)[
+                        0
+                    ][0]
+                except IndexError:
+                    delay = 0
+                delays.append(delay)
+        time_of_flight_readout = float(min(delays) / sampling_rate + MINIMUM_TOF)
         time_of_flights[qubit] = time_of_flight_readout
-
     return TimeOfFlightReadoutResults(time_of_flights)
 
 
@@ -117,16 +143,25 @@ def _plot(
     figures = []
     fitting_report = ""
     fig = go.Figure()
-    qubit_data = data[target]
     sampling_rate = data.sampling_rate
-    y = qubit_data.samples
+    y = magnitude(data.data[target])
 
     fig.add_trace(
         go.Scatter(
             x=np.arange(0, len(y)) * sampling_rate + MINIMUM_TOF,
-            y=y,
+            y=data.data[target][:, 0],
+            textposition="bottom center",
+            name="I",
+            showlegend=True,
+            legendgroup="group1",
+        ),
+    )
+    fig.add_trace(
+        go.Scatter(
+            x=np.arange(0, len(y)) * sampling_rate + MINIMUM_TOF,
+            y=data.data[target][:, 1],
             textposition="bottom center",
-            name="Expectation value",
+            name="Q",
             showlegend=True,
             legendgroup="group1",
         ),
@@ -145,7 +180,19 @@ def _plot(
             line_color="grey",
         )
         fitting_report = table_html(
-            table_dict(target, "Time of flights [ns]", fit.time_of_flights[target])
+            table_dict(
+                target,
+                [
+                    "Intermediate Frequency [Hz]",
+                    "Readout amplitude [a.u.]",
+                    "Time of flights [ns]",
+                ],
+                [
+                    data.intermediate_frequency,
+                    data.amplitude[target],
+                    fit.time_of_flights[target],
+                ],
+            )
         )
     fig.update_layout(
         showlegend=True,

doc/source/protocols/signal/time_of_flight_2d_nodem.png

@@ -0,0 +1,12 @@
+from ...auto.operation import QubitId
+from ...calibration import CalibrationPlatform
+
+
+def _get_lo_frequency(platform: CalibrationPlatform, qubit: QubitId) -> float:
+    """Get LO frequency given QubitId.
+
+    Currently it assumes that instruments with LOs is first one.
+    """
+    probe = platform.channels[platform.qubits[qubit].probe]
+    lo_config = platform.config(probe.lo)
+    return lo_config.frequency