Merge pull request #705 from qiboteam/fix_fit_error

Adding `try-except` block in flux fits

Author: andrea-pasquale

src/qibocal/auto/task.py

@@ -10,7 +10,7 @@
 from qibolab.qubits import QubitId, QubitPairId
 from qibolab.serialize import dump_platform
 
-from ..config import raise_error
+from ..config import log, raise_error
 from ..protocols.characterization import Operation
 from ..utils import (
     allocate_qubits_pairs,
@@ -239,7 +239,12 @@ def update_platform(self, platform: Platform, update: bool):
         """Perform update on platform' parameters by looping over qubits or pairs."""
         if self.task.update and update:
             for qubit in self.task.qubits:
-                self.task.operation.update(self.results, platform, qubit)
+                try:
+                    self.task.operation.update(self.results, platform, qubit)
+                except KeyError:
+                    log.warning(
+                        f"Skipping update of qubit {qubit} due to error in fit."
+                    )
             (self.datapath / PLATFORM_DIR).mkdir(parents=True, exist_ok=True)
             dump_platform(platform, self.datapath / PLATFORM_DIR)
 

src/qibocal/protocols/characterization/flux_dependence/qubit_flux_dependence.py

@@ -12,6 +12,7 @@
 
 from qibocal import update
 from qibocal.auto.operation import Data, Parameters, Qubits, Results, Routine
+from qibocal.config import log
 from qibocal.protocols.characterization.qubit_spectroscopy_ef import (
     DEFAULT_ANHARMONICITY,
 )
@@ -49,8 +50,8 @@ class QubitFluxResults(Results):
     """Sweetspot for each qubit."""
     frequency: dict[QubitId, float]
     """Drive frequency for each qubit."""
-    d: dict[QubitId, float]
-    """Asymmetry."""
+    asymmetry: dict[QubitId, float]
+    """Asymmetry between junctions."""
     fitted_parameters: dict[QubitId, dict[str, float]]
     """Raw fitting output."""
     matrix_element: dict[QubitId, float]
@@ -179,7 +180,7 @@ def _fit(data: QubitFluxData) -> QubitFluxResults:
     qubits = data.qubits
     frequency = {}
     sweetspot = {}
-    d = {}
+    asymmetry = {}
     matrix_element = {}
     fitted_parameters = {}
 
@@ -203,25 +204,33 @@ def _fit(data: QubitFluxData) -> QubitFluxResults:
                 signal,
             )
 
-        popt = curve_fit(
-            utils.transmon_frequency,
-            biases,
-            frequencies * HZ_TO_GHZ,
-            bounds=utils.qubit_flux_dependence_fit_bounds(
-                data.qubit_frequency[qubit], qubit_data.bias
-            ),
-            maxfev=100000,
-        )[0]
-        fitted_parameters[qubit] = popt.tolist()
-        frequency[qubit] = popt[0] * GHZ_TO_HZ
-        d[qubit] = popt[1]
-        sweetspot[qubit] = popt[3]
-        matrix_element[qubit] = popt[2]
+        try:
+            popt = curve_fit(
+                utils.transmon_frequency,
+                biases,
+                frequencies * HZ_TO_GHZ,
+                bounds=utils.qubit_flux_dependence_fit_bounds(
+                    data.qubit_frequency[qubit], qubit_data.bias
+                ),
+                maxfev=100000,
+            )[0]
+            fitted_parameters[qubit] = popt.tolist()
+            frequency[qubit] = popt[0] * GHZ_TO_HZ
+            asymmetry[qubit] = popt[1]
+            sweetspot[qubit] = popt[3]
+            matrix_element[qubit] = popt[2]
+        except ValueError as e:
+            log.error(
+                f"Error in qubit_flux protocol fit: {e} "
+                "The threshold for the SNR mask is probably too high. "
+                "Lowering the value of `threshold` in `extract_*_feature`"
+                "should fix the problem."
+            )
 
     return QubitFluxResults(
         frequency=frequency,
         sweetspot=sweetspot,
-        d=d,
+        asymmetry=asymmetry,
         matrix_element=matrix_element,
         fitted_parameters=fitted_parameters,
     )
@@ -245,7 +254,7 @@ def _plot(data: QubitFluxData, fit: QubitFluxResults, qubit):
                 [
                     np.round(fit.sweetspot[qubit], 4),
                     np.round(fit.frequency[qubit], 4),
-                    np.round(fit.d[qubit], 4),
+                    np.round(fit.asymmetry[qubit], 4),
                     np.round(fit.matrix_element[qubit], 4),
                 ],
             )
@@ -257,7 +266,7 @@ def _plot(data: QubitFluxData, fit: QubitFluxResults, qubit):
 def _update(results: QubitFluxResults, platform: Platform, qubit: QubitId):
     update.drive_frequency(results.frequency[qubit], platform, qubit)
     update.sweetspot(results.sweetspot[qubit], platform, qubit)
-    update.asymmetry(results.d[qubit], platform, qubit)
+    update.asymmetry(results.asymmetry[qubit], platform, qubit)
 
 
 qubit_flux = Routine(_acquisition, _fit, _plot, _update)

src/qibocal/protocols/characterization/flux_dependence/resonator_flux_dependence.py

@@ -11,6 +11,7 @@
 
 from qibocal import update
 from qibocal.auto.operation import Data, Parameters, Qubits, Results, Routine
+from qibocal.config import log
 
 from ..utils import GHZ_TO_HZ, HZ_TO_GHZ, table_dict, table_html
 from . import utils
@@ -38,16 +39,16 @@ class ResonatorFluxResults(Results):
     """Readout frequency for each qubit."""
     sweetspot: dict[QubitId, float]
     """Sweetspot for each qubit."""
-    d: dict[QubitId, float]
-    """Asymmetry."""
+    asymmetry: dict[QubitId, float]
+    """Asymmetry between junctions."""
     bare_frequency: dict[QubitId, float]
     """Resonator bare frequency."""
     drive_frequency: dict[QubitId, float]
     """Qubit frequency at sweetspot."""
     fitted_parameters: dict[QubitId, dict[str, float]]
     """Raw fitting output."""
-    g: dict[QubitId, float]
-    """Coupling."""
+    coupling: dict[QubitId, float]
+    """Qubit-resonator coupling."""
     matrix_element: dict[QubitId, float]
     """C_ii coefficient."""
 
@@ -169,12 +170,12 @@ def _fit(data: ResonatorFluxData) -> ResonatorFluxResults:
     qubits = data.qubits
     frequency = {}
     sweetspot = {}
-    d = {}
+    asymmetry = {}
     bare_frequency = {}
     drive_frequency = {}
     fitted_parameters = {}
     matrix_element = {}
-    g = {}
+    coupling = {}
 
     for qubit in qubits:
         qubit_data = data[qubit]
@@ -196,36 +197,46 @@ def _fit(data: ResonatorFluxData) -> ResonatorFluxResults:
                 signal,
             )
 
-        popt = curve_fit(
-            utils.transmon_readout_frequency,
-            biases,
-            frequencies * HZ_TO_GHZ,
-            bounds=utils.resonator_flux_dependence_fit_bounds(
-                data.qubit_frequency[qubit],
-                qubit_data.bias,
-                data.bare_resonator_frequency[qubit],
-            ),
-            maxfev=100000,
-        )[0]
-        fitted_parameters[qubit] = popt.tolist()
-
-        # frequency corresponds to transmon readout frequency
-        # at the sweetspot popt[3]
-        frequency[qubit] = utils.transmon_readout_frequency(popt[3], *popt) * GHZ_TO_HZ
-        sweetspot[qubit] = popt[3]
-        d[qubit] = popt[1]
-        bare_frequency[qubit] = popt[4] * GHZ_TO_HZ
-        drive_frequency[qubit] = popt[0] * GHZ_TO_HZ
-        g[qubit] = popt[5]
-        matrix_element[qubit] = popt[2]
+        try:
+            popt = curve_fit(
+                utils.transmon_readout_frequency,
+                biases,
+                frequencies * HZ_TO_GHZ,
+                bounds=utils.resonator_flux_dependence_fit_bounds(
+                    data.qubit_frequency[qubit],
+                    qubit_data.bias,
+                    data.bare_resonator_frequency[qubit],
+                ),
+                maxfev=100000,
+            )[0]
+            fitted_parameters[qubit] = popt.tolist()
+
+            # frequency corresponds to transmon readout frequency
+            # at the sweetspot popt[3]
+            frequency[qubit] = (
+                utils.transmon_readout_frequency(popt[3], *popt) * GHZ_TO_HZ
+            )
+            sweetspot[qubit] = popt[3]
+            asymmetry[qubit] = popt[1]
+            bare_frequency[qubit] = popt[4] * GHZ_TO_HZ
+            drive_frequency[qubit] = popt[0] * GHZ_TO_HZ
+            coupling[qubit] = popt[5]
+            matrix_element[qubit] = popt[2]
+        except ValueError as e:
+            log.error(
+                f"Error in resonator_flux protocol fit: {e} "
+                "The threshold for the SNR mask is probably too high. "
+                "Lowering the value of `threshold` in `extract_*_feature`"
+                "should fix the problem."
+            )
 
     return ResonatorFluxResults(
         frequency=frequency,
         sweetspot=sweetspot,
-        d=d,
+        asymmetry=asymmetry,
         bare_frequency=bare_frequency,
         drive_frequency=drive_frequency,
-        g=g,
+        coupling=coupling,
         matrix_element=matrix_element,
         fitted_parameters=fitted_parameters,
     )
@@ -254,8 +265,8 @@ def _plot(data: ResonatorFluxData, fit: ResonatorFluxResults, qubit):
                     np.round(fit.bare_frequency[qubit], 4),
                     np.round(fit.frequency[qubit], 4),
                     np.round(fit.drive_frequency[qubit], 4),
-                    np.round(fit.d[qubit], 4),
-                    np.round(fit.g[qubit], 4),
+                    np.round(fit.asymmetry[qubit], 4),
+                    np.round(fit.coupling[qubit], 4),
                     np.round(fit.matrix_element[qubit], 4),
                 ],
             )
@@ -268,8 +279,8 @@ def _update(results: ResonatorFluxResults, platform: Platform, qubit: QubitId):
     update.bare_resonator_frequency(results.bare_frequency[qubit], platform, qubit)
     update.readout_frequency(results.frequency[qubit], platform, qubit)
     update.drive_frequency(results.drive_frequency[qubit], platform, qubit)
-    update.asymmetry(results.d[qubit], platform, qubit)
-    update.coupling(results.g[qubit], platform, qubit)
+    update.asymmetry(results.asymmetry[qubit], platform, qubit)
+    update.coupling(results.coupling[qubit], platform, qubit)
 
 
 resonator_flux = Routine(_acquisition, _fit, _plot, _update)

src/qibocal/protocols/characterization/flux_dependence/utils.py

@@ -52,7 +52,11 @@ def flux_dependence_plot(data, fit, qubit, fit_function=None):
     )
 
     # TODO: This fit is for frequency, can it be reused here, do we even want the fit ?
-    if fit is not None and not data.__class__.__name__ == "CouplerSpectroscopyData":
+    if (
+        fit is not None
+        and not data.__class__.__name__ == "CouplerSpectroscopyData"
+        and qubit in fit.fitted_parameters
+    ):
         params = fit.fitted_parameters[qubit]
         bias = np.unique(qubit_data.bias)
         fig.add_trace(