code: def measure_wf_coefficients(prep_program, coeff_list, reference_state,
                              quantum_resource, variance_bound=1.0E-6):
      
      num_qubits = len(prep_program.get_qubits())
      normalizer_ops = projector_generator(reference_state, reference_state)
      c0_coeff, _, _ = estimate_locally_commuting_operator(
          prep_program, normalizer_ops, variance_bound=variance_bound,
          quantum_resource=quantum_resource)
      c0_coeff = np.sqrt(c0_coeff)
  
      amplitudes = []
      for ii in coeff_list:
          if ii == reference_state:
              amplitudes.append(c0_coeff)
          else:
              bra = list(map(int, np.binary_repr(ii, width=num_qubits)))
              c_ii_op = projector_generator(reference_state, bra)
  
              result = estimate_locally_commuting_operator(
                  prep_program, c_ii_op, variance_bound=variance_bound,
                  quantum_resource=quantum_resource)
              amplitudes.append(result[0] / c0_coeff)
  
      return amplitudes
docs: Measure a set of coefficients with a phase relative to the reference_state
  
      :param prep_program: pyQuil program to prepare the state
      :param coeff_list: list of integers labeling amplitudes to measure
      :param reference_state: Integer of the computational basis state to use as
                              a reference
      :param quantum_resource: An instance of a quantum abstract machine
      :param variance_bound: Default 1.0E-6.  variance of the monte carlo
                             estimator for the non-hermitian operator
      :return: returns a list of reference_state amplitude + coeff_list amplitudes