def METHOD_NAME():
    n_qubits = 4
    qubits = cirq.LineQubit.range(n_qubits)
    circuit_clifford = cirq.Circuit(cirq.X(qubits[3]))
    circuit_general = cirq.Circuit(cirq.CCX(qubits[0], qubits[1], qubits[2]), circuit_clifford)
    def _run_dfe(circuit):
        class NoiseOnLastQubitOnly(cirq.NoiseModel):
            def __init__(self):
                self.qubit_noise_gate = cirq.amplitude_damp(1.0)
            def noisy_moment(self, moment, system_qubits):
                return [
                    moment,
                    cirq.Moment(
                        [
                            self.qubit_noise_gate(q).with_tags(cirq.ops.VirtualTag())
                            for q in system_qubits[-1:]
                        ]
                    ),
                ]
        noise = NoiseOnLastQubitOnly()
        noisy_simulator = cirq.DensityMatrixSimulator(noise=noise)
        _, intermediate_results = dfe.direct_fidelity_estimation(
            circuit, qubits, noisy_simulator, n_measured_operators=None, samples_per_term=1
        )
        return intermediate_results.pauli_traces, intermediate_results.clifford_tableau is not None
    # Run both algos
    pauli_traces_clifford, clifford_is_clifford = _run_dfe(circuit_clifford)
    pauli_traces_general, general_is_clifford = _run_dfe(circuit_general)
    assert clifford_is_clifford
    assert not general_is_clifford
    assert len(pauli_traces_clifford) == 2**n_qubits
    for pauli_trace_clifford in pauli_traces_clifford:
        pauli_trace_general = [x for x in pauli_traces_general if x.P_i == pauli_trace_clifford.P_i]
        assert len(pauli_trace_general) == 1
        pauli_trace_general = pauli_trace_general[0]
        # The code itself checks that the rho_i is either +1 or -1, so here we
        # simply test that the sign is the same.
        assert np.isclose(pauli_trace_general.rho_i, pauli_trace_clifford.rho_i, atol=0.01)